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a  Hurricane  Storm  Surge  Model  (TR  DRP-95-4) 


ISSUE:  Dredged  material  disposal  sites  located 
in  open  water  are  classified  as  either  dispersive 
or  nondispersive  depending  on  whether  local 
water  velocities  are  strong  enough  to  erode  and 
trmsport  dredged  material  from  the  deposited 
mound.  The  Corps  needs  the  capability  to 
predict  stability  of  the  mound  and  long-term 
migration  patterns  of  eroded  material  to 
(1)  identify  acceptable  disposal-site  locations, 
and  (2)  provide  a  quantitative  approach  for 
gaining  site-designation  approval. 

RESEARCH:  The  overall  work-unit  objective 
is  development  of  a  systematic  approach  for 
predicting  the  dispersion  characteristics  of  a 
specific  open-water  disposal  site.  This  objec¬ 
tive  includes  the  following  goals: 

•  Identify  realistic  wind-,  wave-,  tide-, 
and  storm-generated  velocity  boundary 
conditions. 

•  Develop  numerical  models  capable  of 
simulating  dispersion  characteristics  of 
dredged-material  mounds  for  periods  of 
time  in  excess  of  one  year. 

•  Provide  site-designation  technology  to 
field  engineers  as  a  tool  in  site  identifica¬ 
tion  and  designation. 


One  of  the  areas  of  interest  concerns  the  effect 
of  tropical  storm  surges  on  disposal  site  stability. 
This  report  presents  results  of  studies  undertaken 
to  investigate  the  influence  of  numerical  model 
grid  domain  size  and  grid  stmcture  on  the  accu¬ 
racy  of  numerical  simulations  of  hurricane  surge 
in  the  coastal  region. 

SUMMARY:  The  numerical  model  ADCIRC- 
2DDI  was  used  to  develop  a  relationship  among 
computational  domain  size,  boundaiy  condition 
specification,  and  the  resulting  physics  asso¬ 
ciated  with  hurricane  storm  surge  generation. 
Results  of  the  study  are  guidelines  which 
can  be  used  to  select  a  numerical  domain 
and  subsequently  construct  an  appropriate 
computational  grid. 

AVAILABILITY  OF  REPORT:  The  report  is 
available  through  the  Interlibrary  Loan  Service 
from  the  U.S.  Army  Engineer  Waterways  Exper¬ 
iment  Station  (WES)  Library,  telephone  (601) 
634-2355.  National  Technical  Information  Ser¬ 
vice  (NTIS)  numbers  may  be  requested  from 
WES  librarians. 
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Summary 


This  report  details  two  convergence  studies  undertaken  to  investigate  the 
influence  of  model  domain  size  and  grid  stmcture  on  the  prediction  of  hurri¬ 
cane  storm  surge  in  the  coastal  region.  The  result  of  these  studies  is  a  series 
of  guidelines  which  direct  selection  of  a  domain  size  and  constmction  of  the 
grid  discretization. 

The  domain  size  sensitivity  study  conducted  establishes  a  relationship 
between  domain  size,  boundary  condition  specification,  and  the  resulting  phys¬ 
ics  associated  with  hurricane  storm  surge  generation.  Storm  surge  response 
along  the  Florida  shelf  in  the  Gulf  of  Mexico  due  to  Hurricane  Kate  is  exam¬ 
ined  over  three  domains  using  two  different  open  ocean  boundary  forcing  func¬ 
tions.  The  computed  storm  surge  response  indicates  that  a  small  domain 
situated  primarily  on  the  continental  shelf  is  inadequate  since  cross-shelf 
boundaries  are  in  regions  of  significant  storm  surge  generation  where  surge 
and,  therefore,  boundary  conditions  are  not  known  a  priori.  A  second  domain 
including  the  entire  Gulf  of  Mexico  basin  captures  the  primary  storm  surge 
well  but  may  not  correctly  model  resonant  modes.  The  dependence  of  these 
modes  on  interactions  with  contiguous  basins  makes  accurate  setup  by  the 
boundary  condition  specification  difficult.  The  primary  storm  surge  response 
as  well  as  resonant  modes  excited  by  the  storm  are  best  represented  using  the 
largest  domain  which  encompasses  the  western  North  Atlantic  Ocean,  the 
Caribbean  Sea,  and  the  Gulf  of  Mexico.  This  domain  with  deep  Atlantic 
Ocean  boundaries  facilitates  simple  boundary  condition  specification  and  mini¬ 
mizes  the  influence  of  boundary  conditions  on  storm  surge  generation  in  the 
coastal  region.  Basin  resonant  modes  and  basin-to-basin  interactions  are  also 
captured. 

The  influence  of  grid  structure  on  hurricane  storm  surge  computations  is 
established  by  comparing  storm  surge  elevations  computed  over  14  grid  discre¬ 
tizations  subject  to  4  synthetic  hurricane  forcings.  An  error  analysis  indicates 
that  accurate  predictions  of  storm  surge  result  when  significant  refinement  of 
shoreline  geometry  and  nearshore  regions  is  provided,  along  with  resolution  to 
one  half  the  spatial  scale  of  the  storm  over  deep  waters.  Grid  discretization 
requirements  are  best  met  using  a  graded  grid  structure  which  yields  low  uni¬ 
form  prediction  errors.  The  guidelines  set  forth  for  the  specification  of  a 
domain  size  and  gird  structure  are  applied  relative  to  the  meteorological  forc¬ 
ing  from  Hurricanes  Kate  and  Camille.  The  convergence  analysis  reinforces 


the  necessity  of  using  a  large  domain  which  has  high  levels  of  resolution  in 
near-shore  regions  and  areas  of  complex  coastal  geometry  for  accurate  predic¬ 
tion  of  primary  storm  surge.  In  addition,  easily  implemented  techniques  for 
assessing  grid  performance  relative  to  accurate  storm  surge  prediction  are 
presented. 


1  introduction 


Numerical  modeling  has  become  an  important  tool  for  assessing  the  hydro¬ 
dynamics  of  continental  margin  waters.  It  is  important  to  recognize  that  the 
computed  response  of  these  waters  is  controlled  by  the  various  components 
which  make  up  a  model,  including:  the  governing  equations,  the  boundary 
conditions,  the  forcing  functions,  the  numerics,  grid  structure,  and  the  compu¬ 
tational  domain  itself.  The  more  that  is  understood  about  a  numerical  model 
and  how  its  components  influence  computations,  the  more  successful  the 
model  will  be  in  representing  hydrodynamic  processes  within  shallow  waters. 
The  research  presented  herein  focuses  on  the  numerical  modeling  of  storm 
surge  generation  and  propagation  over  the  continental  shelf  and  in  coastal 
regions  due  to  the  passage  of  a  hurricane.  In  particular,  the  domain  size  and 
grid  structure  associated  with  the  numerical  model  are  examined  in  terms  of 
their  impact  on  computed  storm  surge  response.  Results  from  the  studies 
undertaken  provide  guidelines  for  selection  of  a  domain  size  and  construction 
of  a  grid  discretization  that  are  applicable  to  the  numerical  prediction  of  storm 
surge. 

Storm  surge  is  a  long-period  wave  caused  by  extreme  wind  and  pressure 
forcing.  The  period,  wavelength,  and  amplitude  characteristics  of  the  storm 
surge  depend  on  geometric  properties  of  the  water  body  as  well  as  characteris¬ 
tics  of  the  meteorological  forcing  (Dendrou,  Moore,  and  Myers  1985).  Water 
heights  associated  with  the  storm  surge  are  superimposed  on  water  levels  gen¬ 
erated  by  tidal  forcings.  At  times  when  peak  tidal  and  surge  elevations  coin¬ 
cide,  water  heights  can  become  extreme,  producing  dangerous  conditions  in 
coastal  areas. 

The  impact  of  hurricane  storm  surge  is  well-documented  throughout  history. 
In  recent  memory,  two  particularly  destructive  storms  occurred  during  the 
hurricane  season  of  1992,  Hurricane  Andrew  in  Florida  and  Hurricane  Iniki  in 
Hawaii.  Apart  from  the  number  of  lives  lost  during  these  hurricanes,  the  costs 
due  to  material  damage  and  disruption  to  local  economies  totalled  into  the 
billions  of  dollars.  The  devastation  wrought  by  these  storms  reinforces  the 
seriousness  of  coastal  flooding  in  seaside  communities  caused  by  hurricane 
storm  surge. 

Information  regarding  probable  storm  surge  elevations  generated  by  hurri¬ 
cane  wind  and  pressure  fields  is  desirable  in  a  variety  of  applications.  The 
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National  Weather  Service  (NWS)  issues  warnings  when  hurricanes  threaten 
U.S.  shorelines.  In  defining  the  hurricane  warning  areas,  the  NWS  considers 
the  range  of  expected  peak  surge  heights  which  have  been  extracted  from 
numerical  storm  surge  model  predictions  (Jelesnianski  1979,  Jelesnianski  and 
Taylor  1973,  Jarvinen  and  Lawrence  1985).  Furthermore,  statistical  databases 
of  storm  surge  response  along  the  eastern  U.S.  and  Gulf  of  Mexico  coasts  are 
created  using  water  elevations  computed  by  numerical  storm  surge  model 
hindcasts  of  historical  hurricanes  (Scheffher  et  al.  1994).  In  another  appli¬ 
cation,  design  criteria  for  offshore  oil  structures  are  formulated  with  the  aid  of 
storm  surge  elevation  and  velocity  fields  generated  by  numerical  model  calcu¬ 
lations.  Storm  surge  models  are  also  used  to  develop  coastal  protection  sys¬ 
tems  and  to  investigate  the  fate  of  disposed  dredged  material  within  the  coastal 
environment.  The  success  of  each  of  these  endeavors  rests  upon  accurate 
prediction  of  the  storm  surge  response  due  to  hurricane  forcing  and  consequen¬ 
tly  on  the  formulation  of  an  accurate  numerical  storm  surge  model. 

Three  important  considerations  in  the  numerical  modeling  of  storm  surge 
are  the  adequacy  of  the  grid  resolution,  the  specification  of  boundary  condi¬ 
tions,  and  the  representation  of  resonant  modes  and  surge  forerunner  (Reid 
1990).  The  primary  objective  of  the  research  presented  is  to  address  these 
issues  in  the  context  of  two  sensitivity  studies  examining  the  influence  of 
domain  size  and  grid  structure  on  the  numerical  computation  of  hurricane 
storm  surge.  Findings  from  this  woik  form  a  set  of  guidelines  which  address 
domain  size  selection  and  grid  construction  for  storm  surge  model  applications. 

As  noted  by  Jarvinen  and  Lawrence  (1985),  the  initial  design  and  location 
of  the  basin  domain  are  critical  steps  in  the  formulation  of  a  storm  surge 
model.  Judicious  location  of  a  domain  minimizes  the  effect  of  boundary  con¬ 
dition  specification  and  allows  the  proper  set  of  resonant  modes.  Storm  surge 
modeling  efforts  to  date  have  primarily  used  domains  which  are  limited  to  the 
continental  shelf  region.  Little  effort  has  been  made  towards  assessing  the 
appropriateness  of  such  a  domain  in  storm  surge  models  or  understanding  the 
effects  of  model  domain  size  on  accurate  storm  surge  prediction.  Through  a 
comprehensive  investigation  of  domain  size  sensitivity  in  storm  surge  models 
to  boundary  condition  specification  and  primary  surge  and  surge  forerunner 
generation,  a  domain  size  which  accurately  captures  the  storm  surge  response 
is  established. 

Since  numerical  models  perform  calculations  at  discrete  points,  the  resolu¬ 
tion  provided  by  the  domain  discretization  significantly  influences  computed 
storm  surge  elevations.  An  examination  of  the  relationships  among  grid  spac¬ 
ing,  coastline  variation  and  resolution,  hurricane  forcing,  and  errors  in  the 
computed  storm  surge  provide  insight  into  the  selection  of  an  optimal  grid 
structure.  An  appropriate  grid  structure  uses  the  least  number  of  discrete 
points  to  produce  storm  surge  predictions  having  minimal,  uniform  errors 
throughout  the  domain.  In  the  process  of  studying  the  influence  of  grid  resolu¬ 
tion,  a  simple  yet  effective  method  of  assessing  numerical  errors  associated 
with  the  grid  discretization  is  implemented. 
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Prior  to  the  study  of  domain  size  and  grid  structure,  a  hydrodynamic  model 
used  for  storm  surge  computations  is  selected  in  Chapter  2.  During  the  course 
of  development,  storm  surges  and  their  associated  forenmners  encompass  a 
large  range  of  spatial  scales.  An  acceptable  numerical  storm  surge  model  must 
represent  the  dissipative  processes  acting  on  all  scales  of  motion,  so  that 
energy  in  the  system  is  not  artificially  excited,  and/or  trapped  within  the  sys¬ 
tem,  and/or  aliased  in  nonphysical  ways.  Furthermore,  an  appropriate  storm 
surge  model  incorporates  complex  coastal  geometries,  accounts  for  rapidly 
varying  bathymetry  in  the  continental  slope  and  shelf  regions,  and  permits 
reasonable  boundary  condition  specification. 

A  finite-element-  (FE)-based  numerical  model  with  its  inherent  grid  flexi¬ 
bility  is  selected  as  the  ideal  formulation  for  a  storm  surge  model.  The 
FE-based  storm  surge  model  ADCIRC-2DDI,  the  depth-integrated  portion  of  a 
system  of  two-  and  three-dimensional  codes  named  ADCIRC  (ADvance 
CIRCulation  model  for  shelves,  coasts,  and  estuaries)  developed  by  Luettich, 
Westerink,  and  Scheffher  (1992)  is  utilized.  ADCIRC-2DD1  implements  the 
generalized  wave-continuity  equation  (GWCE)  and  momentum  balance  equa¬ 
tions.  Accuracy  of  the  equations  used  within  ADCIRC-2DDI  is  well- 
documented  with  respect  to  the  solution  of  various  shallow-water  problems 
(Foreman  1988;  Lynch  et  al.  1988;  Walters  1988;  Werner  and  Lynch  1989; 
Walters  and  Werner  1989;  Gray  1989;  Lynch  and  Werner  1991;  Luettich, 
Westerink,  and  Scheffher  1992^  In  addition,  the  finite  element  fonnulation  of 
ADCIRC-2DDI  leads  to  easy  incorporation  of  coastline  detail  and  allows  nodal 
densities  which  range  over  three  to  four  orders  of  magnitude.  The  wide  varia¬ 
tion  in  nodal  density  permits  significant  resolution  of  shoreline  detail  as  well 
as  high  levels  of  refinement  near  coastal  areas  and  in  regions  of  shallow  bathy¬ 
metry  and  rapid  bathymetric  change.  Moreover,  the  discrete  problem  remains 
well  within  computational  limits  despite  the  large  variation  in  nodal  density. 

Aside  from  the  numerical  storm  surge  model,  careful  consideration  is  given 
in  Chapter  3  to  the  specification  of  wind  stress  and  pressure  forcing  at  the 
water  surface.  Wind  stress  and  pressure  forcing  are  what  drive  the  numerical 
storm  surge  model.  As  such,  storm  surge  predictions  are  often  limited  by  the 
accuracy  of  the  specified  wind  stress  and  pressure  fields  (Hubbert,  Leslie,  and 
Manton  1990;  Dendrou,  Moore,  and  Myers  1985;  Flather  1984).  A  numerical 
model  of  the  dynamics  of  the  planetary  boundary  layer,  which  takes  advantage 
of  available  meteorological  data  from  historical  storms,  is  chosen  to  compute 
the  wind  stress  and  pressure  forcing  (Cardone,  Greenwood,  and  Greenwood 
1992).  This  physically  based  wind  model  offers  a  more  realistic  representation 
of  hurricane  winds  than  more  common  empirical  models. 

Upon  establishing  an  acceptable  framewoik  for  the  numerical  storm  surge 
model,  the  domain  size  sensitivity  study  is  conducted  in  Qiapter  4.  For  this 
study.  Hurricane  Kate  (November  1985),  which  made  landfall  on  the  Florida 
shelf  in  the  Gulf  of  Mexico,  is  applied  over  three  domains.  Each  domain 
extends  out  from  the  landfall  point  of  the  hurricane  and  covers  successively 
larger  regions.  The  first  domain  is  a  relatively  small  coastal  domain  which 
extends  mainly  over  the  continental  shelf  in  regions  of  shallow  water.  The 
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second  domain  covers  the  entire  Gulf  of  Mexico  basin.  The  final  domain 
extends  well  into  the  deep  western  North  Atlantic  Ocean.  Storm  surge  compu¬ 
tations  over  these  three  domains  are  compared  to  determine  the  most  appropri¬ 
ate  domain  size  for  accurate  representation  of  the  stonn  surge  response. 
Verification  of  the  storm  surge  predictions  computed  over  the  domain  selected 
as  most  appropriate  with  respect  to  boundary  condition  specification  and  accu¬ 
rate  generation  of  primary  surge  and  resonant  modes  is  achieved  through  a 
comparison  to  measured  storm  surge  elevations. 

An  examination  of  grid  structure  and  its  influence  on  storm  surge  computa¬ 
tions  in  Chapter  5  proceeds  by  considering  idealized  domains  and  hurricane 
forcings.  The  effects  of  grid  spacing,  coastline  variability,  and  coastline  reso¬ 
lution,  as  well  as  the  role  of  various  hurricane  forcing  parameters  on  storm 
surge  generation  are  determined  separately  using  this  simplified  problem.  A 
discussion  of  errors  in  the  storm  surge  predictions  computed  using  four  hurri¬ 
cane  forcings  over  fourteen  grid  discretizations,  two  coastline  representations, 
and  three  discrete  resolutions  of  the  coastline  highlights  important  aspects  of 
the  grid  discretization  necessary  for  accurate  storm  surge  computations. 
Namely,  representation  of  shoreline  geometry  and  refinement  in  near  coastal 
regions  dramatically  influences  computed  storm  surge  response.  In  addition, 
discretization  requirements  over  the  deep  ocean  vary  significantly  from  the 
degree  of  refinement  needed  in  coastal  areas.  Incorporation  of  observations 
regarding  the  effects  of  grid  resolution  leads  to  the  selection  of  a  graded  grid 
structure  appropriate  for  storm  surge  computations.  Furthermore,  the  error 
analysis  techniques  implemented  in  the  study  of  grid  structure  provide  a  practi¬ 
cal  means  for  evaluating  the  perfomiance  of  a  stomi  surge  model  in  relation  to 
its  grid  discretization. 

In  Chapter  6,  an  application  of  the  guidelines  set  forth  for  selection  of  a 
domain  size  and  grid  structure  are  implemented  in  the  context  of  modeling  the 
storm  surge  generated  by  historical  hurricanes  Kate  (November  1985)  and 
Camille  (August  1969).  Use  of  a  large  domain  having  deep  ocean  boundaries 
in  conjunction  with  a  graded  grid  structure  which  provides  significant  refine¬ 
ment  of  shoreline  detail  and  extensive  resolution  in  coastal  areas  is  validated 
by  these  applications.  Furthermore,  these  applications  establish  a  procedure  by 
which  storm  surge  model  predictions  can  be  evaluated  relative  to  the  grid 
discretization  and  a  basis  for  further  grid  refinement  is  provided. 

A  final  statement  in  Chapter  7  summarizes  the  major  conclusions  of  this 
research  and  reiterates  the  implications  of  this  work  on  the  formulation  of 
numerical  storm  surge  models.  Also  mentioned  are  several  issues  relating  to 
the  numerical  prediction  of  storm  surge  which  are  to  be  considered  in  future 
work. 
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2  Hydrodynamic  Model 


The  purpose  of  a  hydrodynamic  model  is  to  characterize  the  important  flow 
features  of  surface  waters  driven  by  tides,  wind,  and  atmospheric  pressure 
gradients.  Complex  circulation  patterns  on  the  continental  shelf  and  in  coastal 
regions  necessitate  the  use  of  a  numerical  hydrodynamic  model.  The  numeri¬ 
cal  hydrodynamic  model  computes  spatial  and  temporal  distributions  of  veloc¬ 
ity  and  sea  surface  elevation  from  which  circulation  patterns  can  be  inferred. 

The  success  of  the  investigation  into  the  characteristics  of  storm  surge 
response  in  the  coastal  ocean  depends  in  part  on  the  selection  of  an  appropriate 
hydrodynamic  model.  An  optimal  model  is  one  that  captures  the  physics  of 
the  storm  surge  response,  yet  remains  computationally  feasible.  In  this  chap¬ 
ter,  the  criteria  used  in  selecting  a  hydrodynamic  model  are  outlined.  In  par¬ 
ticular,  appropriate  governing  equations  are  set  forth  as  well  as  stipulations 
naming  accuracy,  efficiency,  and  grid  flexibility  as  necessary  characteristics  of 
a  numerical  solution  algorithm.  An  overview  and  evaluation  of  current  model¬ 
ing  strategies  is  presented  in  light  of  these  model  selection  criteria.  Subse¬ 
quently,  a  hydrodynamic  model  is  chosen  and  a  description  of  the  model 
follows.  Highlights  and  important  features  of  the  governing  equations  and  the 
numerical  solution  algorithm  implemented  within  the  selected  hydrodynamic 
model  are  also  presented. 


Model  Selection 

Criteria  for  seiection 

The  theoretical  basis  of  a  hydrodynamic  model  is  found  in  the  principles  of 
mass  and  momentum  conservation.  For  vertically  weU-mixed  surface  waters 
experiencing  tidal  and  atmospheric  forcing,  the  flow  physics  are  described  by 
the  shaUow-water  equations,  a  depth-integrated  form  of  the  conservation  laws 
(Le  Mehaute  1976).  The  shallow-water  equations  have  been  used  successfully 
by  engineers  and  researchers  for  years  to  predict  tidal  and  wind-forced  circula¬ 
tions  (e.g.  Johns  and  Ali  1980;  Gray  and  Kinnmark  1983;  Flather  1984; 
Kowalik  1984;  LeProvost  and  Vincent  1986;  Murty,  Flather,  and  Henry  1986; 
Dube,  Sinba,  and  Roy  1986;  Baptista,  Westerink,  and  Turner  1989;  Westerink, 
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Stolzenbach,  and  Connor  1989;  Westerink  et  al.  1992a;  Chen,  Shaffer,  and 
Kim  1993;  Westerink  et  al.  1994a;  Luettich  and  Westerink,  in  preparation). 

Modeling  storm  surge  response  due  to  hurricane  wind  and  pressure  forcing 
in  the  coastal  ocean  is  the  primary  focus  of  this  work.  In  the  shallow  waters 
of  the  continental  shelf,  the  large  magnitude  of  hurricane  wind  forcing  pro¬ 
duces  a  boundary  layer  at  the  sea  surface  which  extends  essentially  throughout 
the  entire  water  column.  As  a  result,  variations  in  the  velocity  over  the  ocean 
depth  are  assumed  negligible,  making  the  depth-averaged  shallow-water  equa¬ 
tions  an  appropriate  framewoik  for  the  hydrodynamic  model. 

Obviously,  the  hydrodynamic  model  formulation  and  solution  strategy  must 
be  accurate.  Wavelength  and  phase  propagation  characteristics,  mass  conserva¬ 
tion  properties,  and  performance  for  test  cases  and  field  applications  are  all 
used  in  the  evaluation  of  model  accuracy. 

Efficiency  of  the  numerical  solution  algorithm  is  the  second  consideration 
in  model  selection.  The  storm  surge  simulations  conducted  herein  occur  over 
several  domain  sizes  including  a  very  large  domain  which  covers  both  conti¬ 
nental  shelf  and  coastal  regions  as  well  as  portions  of  the  deep  ocean.  The 
time  firame  for  these  simulations  ranges  from  days  to  months,  when  tidal  forc¬ 
ing  is  included.  The  feasibility  of  a  model  which  utilizes  large  domains  and 
extended  simulation  periods  depends  upon  the  efficiency  of  the  numerical 
solution  algorithm.  Efficiency  within  a  numerical  algorithm  is  achieved 
through  a  minimization  of  the  number  of  degrees  of  freedom  and  the  number 
of  operations  per  degree  of  freedom  at  each  time-step. 

Lastly,  both  efficiency  and  accuracy  of  the  hydrodynamic  model  formula¬ 
tion  are  enmeshed  with  the  grid  flexibility  characteristics  inherent  in  the 
numerical  solution  algorithm.  Maximum  grid  flexibility  allows  minimization 
of  the  number  of  degrees  of  freedom  while  simultaneously  providing  the  local¬ 
ized  resolution  needed  for  accuracy  of  model  predictions. 

The  selection  criteria  for  the  hydrodynamic  model  are  summarized  as  fol¬ 
lows:  the  governing  equations  are  to  be  based  on  the  shallow-water  equations 
and  the  numerical  solution  algorithm  must  be  highly  efficient,  accurate,  and 
maintain  a  significant  level  of  grid  flexibility. 


Current  modeling  approaches 

The  search  for  accurate  and  efficient  solution  algorithms  for  the  shallow- 
water  equations  has  led  to  a  variety  of  equation  formulations  and  the  use  of 
several  numerical  solution  techniques.  Two  prevalent  numerical  discretization 
strategies  are  the  finite  difference  and  finite  element  methods.  Finite  differ¬ 
ence  methods  discretize  derivative  operators  within  the  model  equations  using 
point  difference  expressions  (Lapidus  and  Finder  1982).  Finite  element  tech¬ 
niques  approximate  the  solution  to  model  equations  through  the  use  of  interpo¬ 
lation  functions  (Celia  and  Gray  1992).  A  review  of  some  finite  element 
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solution  techniques  for  the  shaUow-water  equations  is  given  by  Lee  and 
Froehlich  (1986).  Westerink  and  Gray  (1991),  in  their  recent  review,  note 
increasing  similarities  between  finite-difference-based  and  finite-element-based 
models,  especially  with  regard  to  the  accurate  propagation  of  short  wave¬ 
lengths,  grid  flexibility,  and  complex  boundary  representation. 

Early  finite  difference  solutions  of  the  shallow-water  equations  were  quite 
accurate  owing  in  part  to  the  implementation  of  a  staggered  grid  approach 
which  successfully  avoids  the  introduction  of  artificial  short  waves  (Westerink 
and  Gray  1991).  In  contrast,  early  finite  element  schemes  were  plagued  by 
spurious  oscillations  or  numerical  noise  due  to  the  introduction  of  2  •  Ax 
wavelengths,  a  consequence  of  folded  dispersion  characteristics  (Lynch  1983). 
Note  that  2  •  Ax  wavelengths  are  the  smallest  wavelengths  captured  by  a  mesh 
of  minimum  spacing  Ax.  Finite  element  algorithms  included  excessive  damp¬ 
ing  mechanisms  to  counter  the  generation  of  spurious  modes  but,  as  a  result, 
produced  inaccurate  solutions  (Gray  1982).  It  wasn’t  until  Lynch  and  Gray 
(1979)  introduced  the  wave  continuity  equation  (WCE)  formulation  for  the 
shallow-water  equations  that  the  viability  of  finite  element  approaches 
improved. 

The  WCE  formulation  simply  involves  rearrangement  of  the  shallow-water 
equations  prior  to  spatial  discretization.  Solutions  using  the  WCE  successfully 
suppress  short  wavelengths  without  resorting  to  nonphysical  dissipation.  An 
understanding  of  the  origin  and  behavior  of  spurious  oscillations  in  finite  ele¬ 
ment  formulations  (i.e.,  Walters  and  Carey  (1983))  coupled  with  extensive 
numerical  testing  of  the  WCE  formulation  (i.e.,  Kinnmark  and  Gray  (1984, 
1985);  Luettich,  Westerink,  and  Scheffher  (1992);  Westerink  et  al.  (1992b); 
Kolar  ( 1994a, b))  has  resulted  in  very  accurate  WCE-based  finite  element  solu¬ 
tions  (i.e..  Lynch  and  Gray  (1979);  Foreman  1983). 

In  addition  to  accuracy,  efficiency  of  a  numerical  algorithm  must  also  be 
considered.  Representation  of  complex  circulation  patterns  by  a  hydrodynamic 
model  often  requires  highly  refined  grids  having  very  small  nodal  spacings. 

As  a  consequence,  the  stability  criteria  for  various  time-stepping  schemes 
dictate  use  of  a  very  small  time  interval  for  computations.  The  practicality  of 
using  small  time-steps  for  long-period  simulations  depends  on  the  efficiency  of 
the  numerical  solution  algorithm.  Implicit  time-stepping  schemes  are  generally 
more  stable  and  allow  larger  time-steps,  but  they  result  in  time-dependent 
matrices  which  require  reassembly  and  re-solution  at  every  time-step.  This 
procedure  is  computationally  intensive.  Many  finite  difference  models  over¬ 
come  this  problem  by  implementing  the  alternating  direction  implicit  (ADI) 
solution  algorithm,  which  reduces  a  two-dimensional  problem  to  a  sequence  of 
one-dimensional  problems,  thus  significantly  reducing  the  computational  effort 
(Leendertse  1987).  The  ADI  approach,  however,  cannot  be  applied  to 
unstructured  grid  algorithms  such  as  the  finite  element  solution  strategies. 

In  contrast  to  the  standard  primitive  shallow-water  equation  finite  element 
formulation,  a  WCE-based  finite  element  fonnulation  leads  to  decoupled  eleva¬ 
tion  and  velocity  solutions  and  sparse,  symmetric  matrices  (Lynch  and  Gray 
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1979).  Solutions  for  elevation  are  time-independent  due  to  a  reformulation  of 
the  WCE  by  Kinnmark  (1984)  into  a  GWCE.  Matrices  resulting  from  an 
application  of  the  finite  element  method  to  the  momentum  equations  remain 
time-dependent.  In  the  solution  procedure,  this  time  dependence  is  overcome 
by  lumping  the  matrices,  which  causes  only  a  slight  degradation  in  accuracy. 
As  a  consequence,  velocities  are  computed  by  solving  a  trivial  diagonal  system 
of  equations  (Lynch  and  Gray  1979). 

One  final  consideration  in  selecting  a  numerical  solution  algorithm  is  grid 
flexibility.  The  finite  difference  method  is  not  readily  amenable  to  providing 
high  levels  of  refinement  in  localized  areas.  The  restrictions  on  acceptable 
grid  skewness  and  maximum  cell-to-ceU  size  ratios  (i.e..  Heath,  Johnson,  and 
Kim  (1990);  Celia  and  Gray  1992)  often  result  in  finite  difference  grids  which 
are  overrefined,  causing  an  unnecessary  computational  burden.  To  represent 
detailed  coastline  features,  coordinate  transformation  techniques  have  been 
developed  in  the  context  of  finite  difference  schemes,  but  these  methods  still 
do  not  allow  for  increased  resolution  in  localized  areas. 

In  contrast,  the  advantage  of  the  finite  element  method  lies  in  its  tremen¬ 
dous  grid  flexibility  in  representing  the  complexity  of  coastline  detail  and 
providing  varying  degrees  of  resolution  throughout  the  model  domain.  In 
particular,  finite  element  schemes  based  on  triangular  elements  provide  optimal 
flexibility  in  achieving  local  refinement  (Luettich,  Westerink,  and  Scheffner 
1992).  With  the  ever-increasing  size  and  complexity  of  shallow-water  prob¬ 
lems,  the  degree  of  grid  flexibility  in  the  discretization  strategy  directly  affects 
the  efficiency  and  accuracy  of  the  hydrodynamic  model. 


Selected  hydrodynamic  model 

A  careful  evaluation  of  the  accuracy,  efficiency,  and  grid  flexibility  charac¬ 
teristics  of  current  numerical  solution  algorithms  for  the  shallow-water 
equations  in  the  previous  section  has  led  to  selection  of  a  GWCE  finite- 
element-based  formulation  of  the  shallow-water  equations.  One  hydrodynamic 
model  implementing  a  GWCE  finite-element  formulation  is  ADCIRC-2DDI, 
the  depth-integrated  portion  of  a  system  of  two-  and  three-dimensional  codes 
developed  by  Luettich,  Westerink,  and  Scheffner  (1992)  and  Westerink  et  al. 
(1992b).  All  tidal  and  storm  surge  simulations  conducted  within  utilize  the 
ADCIRC-2DDI  model. 

The  efficiency  and  accuracy  of  ADCIRC-2DDI  are  well  understood  due  to 
extensive  numerical  testing  and  analysis  of  the  model  code.  The  algorithms 
that  comprise  ADCIRC-2DDI  effectively  minimize  the  required  number  of 
degrees  of  freedom  for  a  desired  level  of  accuracy,  show  good  stability  charac¬ 
teristics,  generate  no  spurious  artificial  modes,  have  a  minimum  inherent  artifi¬ 
cial  numerical  damping,  efficiently  separate  the  partial  differential  equations 
into  systems  of  algebraic  equations  with  time-independent  matrices,  and  are 
capable  of  running  months  to  years  of  simulation  while  providing  detailed 
computations  of  the  circulation  patterns  within  a  water  body  (Luettich, 
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Westerink,  and  Scheffher  1992).  An  investigation  into  the  mass  conservation 
properties  of  the  GWCE  formulation  by  Kolar  et  al.  (1994b)  provides  an 
understanding  of  the  mass  balance  characteristics  of  the  ADC1RC-2DDI  model. 
In  addition,  Westerink  et  al.  (1994a)  confirm  the  accuracy  of  the  boundary 
condition  formulation  within  ADCIRC-2DDI. 

The  ADCIRC-2DDI  hydrodynamic  model  produces  computations  which  are 
in  agreement  with  two  established  benchmark  problems,  the  quarter  annular 
test  case  (Gray  1982)  and  the  North  Sea/English  Channel  system  (Gray  1989). 
The  ADCIRC-2DDI  model  has  also  been  applied  to  a  number  of  field  studies 
with  excellent  results  (Westerink,  Stolzenbach,  and  Connor  1989;  Westerink 
et  al.  1992a,b;  Westerink,  Luettich,  and  Muccino  1994;  Kolar  et  al.  1994a,b). 
Detailed  documentation  (Luettich,  Westerink,  and  Scheffner  1992)  as  well  as 
an  instructive  commentary  regarding  appropriate  modeling  strategies  for  use 
with  ADCIRC-2DDI  (Westerink  et  al.  1992b)  facilitates  application  of  the 
ADCIRC-2DDI  model. 


ADCIRC-2DDI  Model  Description 

Governing  equations 

The  generalized  wave  continuity  equation  formulation  in  ADCIRC-2DDI  is 
based  on  the  well-known  shallow-water  equations  (Le  Mehaute  1976; 

Kinnmark  1984).  The  primitive  form  of  the  shallow-water  equations  is  derived 
by  averaging  the  conservation  laws  of  mass  and  momentum  over  the  time  scale 
of  turbulent  fluctuations  and  the  ocean  depth.  Turbulent  fluctuations  at  the 
microscale  present  in  all  hydrodynamic  flows  are  characterized  by  spatial  and 
temporal  variations  in  the  velocity  and  pressure  fields  and  are  a  mechanism  for 
momentum  transfer.  An  application  of  time  averaging  to  the  conservation  laws 
provides  a  means  to  parameterize  these  microscale  fluctuations  in  terms  of 
macroscale  quantities.  The  scale  for  time  averaging  is  selected  to  be  long 
enough  to  capture  a  statistically  significant  sample  of  the  turbulent  fluctuations 
and  yet  short  enough  so  that  macroscopic  variations  of  the  average  quantity  are 
not  included  (Gray  et  al.  1993).  Within  the  shallow-water  equations,  turbulent 
fluctuations  are  represented  in  an  average  way  through  the  momentum  diffu¬ 
sion  terms. 

Depth  averaging  of  the  conservation  laws  reduces  a  three-dimensional  prob¬ 
lem  to  a  two-dimensional  one,  leaving  as  unknowns  the  surface  water  elevation 
^  and  the  depth-averaged  lateral  velocities  U  and  V.  For  flows  having  small 
vertical  velocity  gradients  (e.g.  well-mixed  systems)  or  environments  where 
lateral  flows  are  quite  large  in  comparison  to  vertical  velocities  (e.g.  nearly 
horizontal  flow),  an  application  of  the  shallow-water  equations  is  accepted 
practice. 

Derivation  of  the  shallow-water  equations  from  the  time  and  depth-averaged 
conservation  laws  involves  assumptions  of  incompressibility  and  hydrostatic 
conditions  as  well  as  the  Boussinesq  approximation.  In  practice. 
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incompressibility  of  a  fluid  implies  that  density  variations  seen  while  moving 

with  the  fluid  are  negligible  (i.e.,  =  0).  Furthermore,  the  Boussinesq 

approximation  specifies  density  as  a  constant  value  (p  =  p^,)  except  when 
gradients  of  density  are  considered.  In  flows  appropriate  for  two-dimensional 
depth-averaged  equations,  the  vertical  acceleration  is  assumed  negligible. 
Consequently,  the  momentum  equation  over  the  vertical  reduces  to  a  balance 
between  pressure  and  gravitational  forces.  This  is  commonly  referred  to  as  the 
hydrostatic  assumption.  Other  simplifications  implicit  within  the  shallow-water 
equations  include  the  neglect  of  changes  in  the  position  of  the  ocean  floor  with 
respect  to  time  and  the  exclusion  of  mass  exchanges  with  the  environment 
(i.e.,  evaporation,  precipitation,  overland  flow,  and  groundwater  interactions). 

For  a  Cartesian  coordinate  system,  conservative  forms  of  the  shaUow-water 
equations  are  written  as  follows: 

a?  ,  dvH  ^  aw  ,  0  (1) 

dt  dx  dy 


dUH  ^  dUUH  ^  dUVH  ^ 

dt  dx  dy  dx 


+  g(C-ttTl) 


(2) 


+  M  +D  +  B  +  —  -  — 
"  "  Po  Pa 


dVH  ,  dVUH  ^  dWH  ^  d 

dt  dx  dy  dy 


Po 


^(C-ari) 


(3) 


*  M*  D  *  B  ^  ^  -  h. 

’  ’  ’  P,  P. 


where 

^  =  free  surface  elevation  relative  to  the  geoid 
t  =  time 

U  =  depth-averaged  horizontal  velocity 
H  =  water  column  depth,  +  h 
h  =  bathymetric  depth 
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X  =  Cartesian  coordinate  direction 
V  =  depth-averaged  horizontal  velocity 
y  =  Cartesian  coordinate  direction 
/=  Coriolis  parameter 
=  atmospheric  pressure  at  the  free  surface 
Pg  =  reference  density  of  water 
g  =  acceleration  due  to  gravity 
a  =  earth  elasticity  factor 
T|  =  Newtonian  equilibrium  tide  potential 
=  depth-integrated  horizontal  momentum  diffusion 
Djj  =  depth-integrated  horizontal  momentum  dispersion  term 
=  depth-integrated  baroclinic  forcing 
=  applied  free  surface  stress 
'^bx  -  applied  bottom  stress 

My  =  depth-integrated  horizontal  momentum  diffusion 
Dy  =  depth-integrated  horizontal  momentum  dispersion  term 
By  =  depth-integrated  baroclinic  forcing 
=  applied  free  surface  stress 
\  -  applied  bottom  stress 

Further  justification  regarding  the  appropriateness  of  these  equations  in 
modeling  the  flow  of  tidal  and  atmospheric  forces  is  provided  by  Blumberg 
and  Mellor  (1987);  Westerink,  Stolzenbach,  and  Connor  (1989);  and  Luettich 
etal.  (1992). 

In  operator  notation,  a  relationship  between  the  conservative  and  nonconser¬ 
vative  momentum  equations  is  given  (Kolar  et  al.  1994b): 
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(4) 


H 


where 


=  nonconservative  momentum  equations 

=  conservative  momentum  equations  (Equations  2  and  3) 

V  =  horizontal  velocity  vector 

C  =  primitive  continuity  equation  (Equation  1) 

Substitution  of  Equations  1-3  into  Equation  4  leads  to  a  reformulation  of 
the  momentum  equations  into  a  primitive,  nonconservative  form: 


^  -fV  =  aq) 

dt  dx  dy  dx  pg 


T  T. 
sx  bx 


T  r  A  X  X  ^  _ 

^  Po  9o 


^fU  =  ^  g(C-aq) 

dt  dx  dy  dy  p^ 


hJ-M+D  +  B  +  ^  ^ 

^  ^  Po  Po 


A  rigorous  derivation  of  Equations  1-6  is  presented  by  Kolar  et  al.  (1994b) 
and  will  not  be  repeated  here. 

The  implementation  of  a  standard  quadratic  parameterization  for  bottom 
stress  and  the  neglect  of  baroclinic  and  lateral  diffusion/  dispersion  terms  leads 
to  a  modified  form  of  the  primitive,  nonconservative  shallow-water  equations: 

(7) 

dt  dx  dy 
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(8) 


dU 

"ST 


Ps 

Po 


+  g(C-an) 


+  _!L  -  xM 
PoH 


dv 

W 


jj'dV  ,,dV  .jj  d 

+  U—  +  V  +  /t/  = 
dx  ay  ay 


fi 

Po 


+  g(C  -  ari) 


+ 


(9) 


with  X*  given  by  the  expression  Cf(L/^  +  V^j  /  H  for  Cy  equal  to  the  bottom 
friction  coefficient. 

Recall  that  one  of  the  objectives  of  the  ADCIRC-2DD1  hydrodynamic 
model  is  to  consider  very  large  domain  problems.  The  use  of  large  domains 
coupled  with  the  recognition  that  tidal  forcings  are  a  global  phenomenon 
necessitates  the  inclusion  of  effects  caused  by  the  curvature  of  the  Earth's 
surface.  Thus,  the  shallow-water  equations  (Equations  7-9)  are  recast  into 
spherical  coordinates  (Flather  1988,  Kolar  et  al.  1994a): 


aC  ^  1  dUH 

dt  /?COS(j)  dX 


^  3(F//cos(t)) 
3<!) 


=  0 


(10) 


du  ^  u  du  ^  V  du 
~W  Rcoscj)  dX  R  04) 


tdsv^U 

R 


1  0  Ps 

/?COS4)  0X.  Po 


+  g(^  -  Tl) 


(11) 


0V  ^  U  dv  ^  V  dv 
dt  Rcos(p  dX  R  04) 


tan^U 

~R~ 


j_  0  ^ 

Po 


+  ^  (C  -'H) 


Po^ 


(12) 
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where 


R  =  radius  of  the  Earth 

A,,  (|)  =  degrees  longitude  (east  of  Greenwich  positive)  and  degrees 
latitude  (north  of  the  equator  positive) 

/  =  2Qsin^ 

=  angular  speed  of  the  Earth 
=  applied  free  surface  stresses 

A  practical  expression  for  the  Newtonian  equilibrium  tide  potential  is  given 
by  Reid  (1990)  as: 


niXM  =  Yf  CjJj„(tQ)Lj(<^)cos[2K(t  -  tQ)/Tj„ 


(13) 


where 


t  =  time  relative  to  tQ 
Iq  =  reference  time 

Cj„  =  constant  characterizing  amplitude  of  tidal  constituent  n  of  species  j 
fj^  =  time-dependent  nodal  factor 


Tj„  =  period  of  constituent  n  of  species  j 

Vj^  =  time-dependent  astronomical  argiunent 

j  =  0,  1,  2  =  tidal  species  (j  =  0,  declinational;  j  =  I,  diurnal; 
j  =  2,  semidiurnal 

Values  for  Cj^  are  presented  by  Reid  (1990).  The  value  for  the  Earth  elas¬ 
ticity  factor  is  often  taken  as  0.69  for  all  tidal  constituents  (Schwiderski  1980, 
Hendershott  1981)  although  its  value  has  been  shown  to  be  slightly 
constituent-dependent,  ranging  between  0.693  and  0.736  (Wahr  1981, 
Woodworth  1990). 

The  Earth's  curvature  must  be  accounted  for  not  only  in  the  governing 
equations  but  also  in  the  finite  element  discretization  (Kolar  et  al.  1994a).  The 
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finite  element  method,  as  briefly  mentioned  in  the  section  titled  “Current  Mod¬ 
eling  Approaches,”  approximates  a  solution  using  interpolating  functions. 

These  interpolating  functions  are  defined  over  elements  and  these  elements  are 
most  often  cast  within  the  framework  of  a  Cartesian  coordinate  system.  In 
order  to  conveniently  implement  the  finite  element  method,  the  governing 
equations  in  spherical  coordinates  are  projected  onto  a  planar  surface  using 
cartographic  projection  techniques.  The  mapping  of  spherical  equations  (Equa¬ 
tions  10-12)  to  a  rectilinear  coordinate  system  is  accomplished  using  a  Carte 
ParaUelogramatique  Projection  (CP)  (Pearson  1990): 

y!  =  /?(X  -  X,g)cos(t)^  (14) 


where  (X^,  %)  is  the  center  point  of  the  projection.  An  application  of  the  CP 
projection  to  Equations  10-12  yields  shallow-water  equations  in  primitive, 
nonconservative  form  expressed  in  the  CP  coordinate  system: 


dt,  ^  cos(t)o  d(UH)  ^  1  a(V7/cos(j))  ^ 

"3?  cos<t)  ay  cos$  ^ 


du  ^  cos^o  ^jdu  ^ydu 
"3F  cos<j)  ^ 


COS<i)o  dPs  .y  S  ^sX  JJ 

— 1"3:7  —  ■  ft)  +  — TT  - 

COS<^  dx  Po  Po^ 


^  cos<j)  R 


=  -  ^  g  (C  +  ft)  +  -4  -  t.v 

Po  Po^ 
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As  discussed  in  the  section  titled  “Current  Modeling  Approaches,”  utilizing 
the  finite  element  method  to  resolve  the  spatial  dependence  of  the  shaUow- 
water  equations  in  their  primitive  form  leads  to  inaccurate  solutions  with 
severe  artificial,  near  2  •  Ax  modes  (Gray  1982;  Lynch  1983;  Westerink, 
Stolzenbach,  and  Connor  1989).  A  reformulation  of  the  primitive  equations 
into  a  GWCE  form  gives  highly  accurate,  noise-free,  finite-element-based 
solutions  to  the  shaUow-water  equations  (Lynch  and  Gray  1979,  Kinnmark 
1984).  The  high  accuracy  of  GWCE-based  FE  solutions  is  a  result  of  their 
excellent  numerical  amplitude  and  phase  propagation  characteristics.  In  fact, 
Fourier  analysis  indicates  that  in  constant-depth  water,  using  linear  interpola¬ 
tion,  a  linear  tidal  wave  resolved  with  25  nodes  per  wavelength  is  more  than 
j^quately  resolved  over  the  range  of  Courant  numbers,  C  =  At  I  Ax,  less 
man  or  equal  to  one  (Luettich,  Westerink,  and  Scheffner  1992).  Furthermore, 
the  monotonic  dispersion  behavior  of  GWCE-based  FE  solutions  avoids  gener¬ 
ating  artificial,  near  2  •  Ax  modes,  which  have  plagued  the  primitive  FE-based 
solutions,  in  both  interior  and  botmdary  equations  (Platzman  1981;  Foreman 
1983;  Westerink  et  al.,  in  preparation  (b)).  Note  that  the  monotonic  dispersion 
behavior  of  GWCE-based  FE  solutions  is  very  similar  to  that  associated  with 
staggered  finite  difference  solutions  to  the  primitive  shaUow-water  equations 
(Westerink  and  Gray  1991).  GWCE-based  FE  solutions  to  the  shallow-water 
equations  allow  for  extremely  flexible  spatial  discretizations  which  result  in  a 
highly  effective  minimization  of  the  discrete  size  of  any  problem  (Le  Provost 
and  Vincent  1986;  Foreman  1988;  Vincent  and  Le  Provost  1988;  Westerink, 
Luettich,  and  Muccino  1994;  Westerink  et  al.,  in  preparation  (a);  Luettich  and 
Westerink,  in  preparation). 

Derivation  of  the  GWCE  is  presented  concisely  by  the  operator  notation 
invoked  by  Kinnmark  (1984)  and  Kolar  (1992): 


GWCE  =  _  +  T  C  -  V  • 
dt 


(19) 


where  GWCE  is  the  generalized  wave  continuity  equation  and  is  a  non¬ 
physical  constant  in  time  and  space  which  controls  the  balance  between  primi¬ 
tive  and  wave  equation  formulations  (Lynch  and  Gray  1979;  Kinnmark  1984; 
Luettich,  Westerink,  and  Scheffner  1992;  Kolar  et  al.  1994b).  Substituting 
continuity  Equation  16  and  the  conservative  forms  of  momentum  Equations  17 
and  18  into  Equation  19,  the  GWCE  in  the  CP  coordinate  system  is: 
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(20) 


The  advective  terms  in  the  GWCE,  written  in  nonconservative  form, 
improve  global  and  local  mass  conservation  (Kolar  et  al.  1994b)  as  well  as 
numerical  stability,  especially  for  advection  dominant  flows  (Westerink  et  al. 
1992b).  The  ADCIRC-2DDI  hydrodynamic  model  solves  the  GWCE,  Equa¬ 
tion  20,  in  conjunction  with  the  primitive  momentum  equations  in  nonconserv¬ 
ative  form.  Equations  17  and  18. 


Numerical  solution 

The  niunerical  discretization  of  the  GWCE  equation  (Equation  20)  and 
momentum  Equations  17  and  18  is  implemented  in  three  stages.  First,  the 
symmetrical  weak  weighted  residual  (SWWR)  statements  for  Equations  17, 

18,  and  20  are  developed.  This  procedure  is  based  on  a  standard  Galerkin 
finite  element  formulation  (Becker,  Carey,  and  Oden  1981;  Celia  and  Gray 
1992).  One  consequence  of  the  SWWR  fonn  is  that  the  order  of  the  deriva¬ 
tives  in  the  governing  equations  is  reduced,  leading  to  a  requirement  of  only 
Cg  functional  continuity  (i.e.,  only  interpolating  functions  themselves,  not  their 
derivatives,  need  be  continuous  between  discrete  points).  Next,  the  stable  and 
accurate  time  discretization  strategies  of  Kinnmark  (1984)  and  Werner  and 
Lynch  (1989)  are  implemented.  A  variably  weighted  three-time-level  implicit 
scheme  is  applied  to  most  linear  terms  in  toe  GWCE.  The  nonlinear  Coriolis, 
atmospheric  pressure,  and  tidal  potential  terms  are  all  treated  explicitly.  Alter¬ 
natively,  convective  terms  within  toe  GWCE  are  evaluated  at  two  known  time 
levels.  This  time  discretization  strategy  results  in  a  linear  algebraic  system  of 
equations  associated  with  toe  GWCE,  which  is  solved  for  unknown  elevations. 
For  toe  momentum  equations,  a  Crank-  Nicolson  two-time-level  implicit 
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scheme  is  used  for  all  terms  except  the  bottom  friction  and  convective  terms, 
which  are  treated  explicitly. 

The  final  step  in  the  numerical  discretization  scheme  is  approximation  of 
the  spatial  domain  using  the  finite  element  method.  Variables  are  expanded 
using  a  Q  interpolation  basis  over  three-node,  linear  triangular  elements. 
Elemental  equations  are  summed  over  the  global  domain  and  inter-element  Q 
functional  continuity  is  enforced.  Details  of  the  finite  element  implementation 
on  a  term-by-term  basis  is  presented  by  Luettich,  Westerink,  and  Scheffner 
(1992)  and  Westerink  et  al.  (1992b).  It  is  noted  that  shoreline  wetting  and/or 
drying  is  not  currently  included  in  the  ADCIRC-2DD1  code. 

The  fully  discretized  model  equations  are  written  in  matrix  notation 
(Luettich,  Westerink,  and  Scheffner  1992): 


y^WCE^k*l  _  pGWCE 

(21) 

•^^XME jj k  +  l  =  p^-^E 

(22) 

j^(t)XW£  ^*  +  1  +  =  p^ME 

(23) 

where 

^wcE  _  5anded,  time-independent  mass  matrix  in 
the  GWCE  equation 

=  vector  of  unknown  surface  elevations  at 
time  level  k  +  \ 

pGwcE  _  vector  of  known  forcings  in  the 
GWCE  equation 

^UME,  -  fime-dependent,  lumped  mass  matrices  in 

the  X,  (j)  directions 

y*  +  1  yE  +  1  _  yectors  of  unknown  velocity  components  at 
the  A:  -I-  1  time  level 

pXME^  p^ME  _  yectors  of  known  forcings  for  the 
momentum  equations 

Elevation  boundary  conditions  are  enforced  within  the  load  vector  of 
the  GWCE  equation  (Equation  21)  and  zero  normal  velocity  boundary  con¬ 
ditions  are  enforced  in  the  momentum  equations  (Equations  22  and  23). 
Westerink  et  al.  (in  preparation  (b))  have  shown  that  solutions  to  the  GWCE 
equation  are  insensitive  to  this  standard  boundary  condition  formulation. 
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The  decoupled  discrete  GWCE  and  momentum  equations  (Equations  21-23) 
lead  to  a  sequential  solution  procedure.  The  GWCE  equation  (Equation  21)  is 
solved  at  each  time-step  for  the  surface  water  elevations,  ^  ^  The  GWCE 
mass  matrix  is  time  independent,  so  it  is  assembled  and  decomposed 

only  once.  The  banded  structure  of  this  matrix  is  not  utilized  by  the  iterative 
preconditioned  conjugate  gradient  matrix  solver  used  for  computations  per¬ 
formed  in  this  work  (Press  et  al.  1986;  Kincaid  and  Cheney  1991).  The  pre¬ 
conditioned  conjugate  gradient  solver  is  selected  because  of  its  efficiency  at 
minimizing  the  memory  requirements  of  ADCIRC-2DDI  for  larger  problems. 
The  load  vector  is  updated  at  each  time-step  with  newly  computed 

surface  water  elevations  and  velocities  from  the  previous  time-step. 

Upon  solving  the  GWCE  for  surface  water  elevations,  the  computed  eleva¬ 
tions  are  substituted  into  momentum  Equations  22  and  23  prior  to  solution  for 
the  velocity  components  '  and  V*  ^  The  time-dependent  mass  matrices 
in  the  momentum  equations  (M^^,  are  lumped  to 

yield  diagonal  matrices  which  require  trivial  solution.  The  lumping  procedure 
applied  here  does  not  introduce  significant  errors  as  shown  by  Lynch  and  Gray 
(1979). 

The  numerical  solution  algorithm  just  described  for  the  ADCIRC-2DDI 
model  is  implemented  in  fuUy  vectorized  form.  A  consequence  of  this  solu¬ 
tion  procedure  is  a  highly  efficient  code  in  terms  of  central  processing  unit 
requirements  per  node.  This  efficiency  is  largely  due  to  the  fact  that  GWCE- 
based  FE  solutions  to  the  shallow-water  equations  allow  for  extremely  flexible 
spatial  discretizations  which  result  in  a  highly  effective  minimization  of  the 
discrete  size  of  any  problem  (Le  Provost  and  Vincent  1986;  Foreman  1988; 
Vincent  and  Le  Provost  1988;  Westerink  et  al.  1992a,b;  Westerink,  Luettich, 
and  Muccino  1994). 


Summary 

Numerical  formulation  of  the  hydrodynamic  model  is  critical  to  the  success 
of  any  storm  surge  modeling  application.  The  depth-averaged  shaUow-water 
equations  are  designated  as  an  appropriate  theoretical  framework  for  the  hydro- 
dynamic  model.  Additional  requirements  of  accuracy,  efficiency,  and  grid 
flexibility  are  imposed  on  the  numerical  solution  technique.  A  brief  review  of 
current  hydrodynamic  modeling  approaches  is  presented  in  terms  of  the  model 
selection  criteria  set  forth.  From  the  assessment  of  available  numerical  solu¬ 
tion  techniques  for  the  shallow- water  equations,  a  GWCE-based  finite  element 
approach  is  chosen  as  the  optimal  solution  strategy.  The  GWCE-based  finite 
element  equations  are  implemented  within  the  hydrodynamic  model 
ADCIRC-2DDI. 

ADCIRC-2DDI  has  demonstrated  the  desired  characteristics  of  accuracy, 
efficiency,  and  a  high  degree  of  grid  flexibility  in  numerous  numerical  tests, 
analyses,  and  field  applications.  Thus,  the  ADCIRC-2DDI  hydrodynamic 
model  is  adopted  for  all  storm  surge  simulations  presented  in  this  work. 
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Details  of  the  governing  equations,  their  limitations,  advantages,  and  numerical 
solution  are  discussed. 

The  ADCIRC-2DDI  hydrodynamic  model  provides  a  sound  tool  with  which 
to  study  the  characteristics  of  the  computed  storm  surge  response  and  to 
address  convergence  issues  related  to  storm  surge  modeling.  When  using  the 
ADCIRC-2DDI  model,  storm  surge  generation  is  simulated  over  large  domains 
for  lengthy  time  periods  with  striking  efficiency  and  accuracy. 
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3  Meteorological  Forcing 


Storm  surge  development  can  be  attributed  to  several  mechanisms 
(Dendrou,  Moore,  and  Myers  1985;  Reid  1990)  each  of  which  is  caused  pri¬ 
marily  by  some  combination  of  wind  and  pressure  forcing.  Furthermore,  it  is 
well-known  that  storm  surge  models  are  highly  sensitive  to  the  representation 
of  wind  forcing  (i.e.,  Chen,  Shaffer,  and  Kim  1993;  Hubbert,  Leslie,  and 
Manton  1990;  Dendrou,  Moore,  and  Myers  1985;  Flather  1984).  Thus,  realis¬ 
tic  storm  surge  generation  by  a  numerical  model  depends  upon  accurate  speci¬ 
fication  of  wind  and  pressure  fields.  Though  measured  values  of  wind  speeds 
and  pressures  can  be  used,  this  data  is  almost  always  scarce  and  too  limited  in 
scope.  More  commonly,  either  an  empirical  or  conservation-based  model  is 
used  to  generate  wind  and  pressure  fields  associated  with  a  hurricane. 

This  chapter  begins  by  offering  a  physical  description  of  the  wind  and 
pressure  forces  within  a  hurricane  and  their  role  in  stonn  surge  generation. 

The  current  state  of  hurricane  models  is  briefly  mentioned  and  selection  of  a 
wind  model  is  made.  The  chosen  wind  model  is  then  described  in  terms  of  its 
governing  equations,  numerical  implementation,  and  interface  with  the  hydro- 
dynamic  model  ADCIRC-2DD1. 

Two  hurricanes  of  record  (Kate  (1985)  and  Camille  (1969))  are  applied  as 
meteorological  forcing  for  the  storm  surge  computations  in  this  work.  For 
each  of  these  hurricanes,  the  path,  wind  speed,  eye  pressure,  and  other  impor¬ 
tant  parameters  are  presented.  A  series  of  synthetic  hurricanes  used  in  the 
investigation  of  grid  structure  influence  on  stonn  surge  generation  are  also 
described. 


Hurricane  Forcing  and  Modeling  Approaches 

A  hurricane  is  a  rotating  cyclone  in  the  tropical  ocean  having  considerable 
size  and  intensity.  More  specifically,  a  hurricane  is  classified  by  wind  speeds 
of  33  m/s  (74  mph)  or  greater  and  is  characterized  by  a  low  pressure  center. 
The  region  of  low  pressure  at  the  eye  of  the  storm  produces  a  rise  in  sea  sur¬ 
face  elevation  named  an  inverted  barometer  effect.  Outwards  from  the  eye, 
pressures  gradually  increase  towards  the  storm's  edge. 
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Hurricane  wind  fields  are  characterized  by  high  winds  blowing  spirally 
around  the  low  pressure  center  of  the  storm.  In  the  Northern  Hemisphere, 
hurricane  wind  fields  orient  in  a  counterclockwise  direction  due  to  the  Coriolis 
forcing,  which  accounts  for  rotation  of  the  Earth.  Relatively  low  wind  speeds 
are  experienced  at  the  hurricane  eye  but  rapidly  increase  out  to  a  radius  rang¬ 
ing  anywhere  from  15  km  to  85  km  (8.1-45.8  N  miles),  a  distance  often 
termed  the  radius  to  maximum  wind.  Wind  speeds  diminish  beyond  the  radius 
to  maximum  wind  out  to  the  storm's  periphery. 

Hurricane  winds  exert  a  shearing  stress  at  the  ocean  surface.  Drag  forces 
associated  with  the  wind  shear  work  both  to  maintain  a  sloping  water  surface 
against  the  shoreline  and  drive  surface  water  in  the  direction  of  the  wind,  the 
latter  being  the  primary  mechanism  for  storm  surge  generation.  Further  details 
regarding  the  formation,  life  cycle,  and  impact  of  a  hurricane  are  given  by 
Dunn  and  Miller  (1964),  Stubbs  (1966),  and  Simpson  and  Riehl  (1981). 

Recent  numerical  models  of  the  hurricane  itself  are  based  on  three- 
dimensional  primitive  conservation  equations  which  are  solved  over  a 
multiply-nested  movable  grid  (Powell  1982;  Bender,  Tuleya,  and  Kurihara 
1984;  Jones  1987;  Kurihara  et  al.  1990).  These  models  are  computationally 
intensive  and  fairly  complex.  Hydrodynamic  modelers  have  traditionally 
sought  simpler  ways  of  estimating  the  wind  and  pressure  forcing  generated  by 
a  hurricane.  In  many  storm  surge  modeling  studies  parameterizations  or 
empirical  models  are  used  to  represent  wind  and  pressure  forcing  (e.g.  Johns 
and  Ali  1980;  Johns  et  al.  1983a,b;  Flather  1984;  Jarvinen  and  Lawrence  1985; 
Shaffer,  Jelenianski,  and  Chen  1986;  Wang  1987;  Westerink,  Stolzenbach,  and 
Connor  1989;  Hearn  and  Holloway  1990).  Murty,  Flather,  and  Henry  (1986) 
review  many  of  the  currently  used  empirical  cyclone  models  and  offer  some 
inter-comparisons.  The  most  complex  empirical  model  reviewed  by  Murty, 
Flather,  and  Henry  (1986)  is  the  Standard  Project  Hurricane  (SPH)  model 
(Cialone  1991)  used  by  Westerink,  Stolzenbach,  and  Connor  (1989),  among 
others.  Westerink,  Stolzenbach,  and  Connor  (1989)  report  severe  overpredic¬ 
tions  of  hurricane  winds  coming  off  the  land  when  using  the  SPH  model.  To 
obtain  satisfactory  results,  an  exponential  decay  factor  was  applied  within  the 
SPH  model  to  winds  computed  near  land. 

To  improve  the  quality  of  the  specified  meteorological  forcing,  a  hurricane 
model  which  accounts  for  the  basic  physical  and  dynamical  characteristics  of  a 
tropical  cyclone  is  sought.  Chen,  Shaffer,  and  Kim  (1993)  suggest  initially 
computing  a  pressure  field  and  then  using  this  pressure  field  to  calculate  geo- 
strophic  winds.  One  model  adopting  such  an  approach  is  the  HURWIN  wind 
model  of  Cardone,  Greenwood,  and  Greenwood  (1992).  The  HURWIN  model 
uses  the  conservation  of  momentum  as  a  foundation  and  is  selected  here  to 
compute  the  meteorological  forcing  for  the  hydrodynamic  model, 
ADCIRC-2DDI.  Though  specification  of  the  wind  and  pressure  forcing  is 
important  to  the  success  of  storm  surge  model  predictions,  any  further  review 
of  hurricane  models  or  testing  of  the  HURWIN  wind  model  is  beyond  the 
scope  of  this  work. 
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HURWIN  Wind  Model 


Governing  equations 

The  HURWIN  wind  model  is  part  of  the  Coastal  Modeling  System  at  the 
U.S.  Army  Corps  of  Engineers  Waterways  Experiment  Station  and  is  based  on 
the  vortex  model  of  Chow  (1971)  modified  by  Cardone,  Greenwood,  and 
Greenwood  (1992).  The  HURWIN  wind  model  is  a  diagnostic  model  and  thus 
is  not  a  predictive  tool.  As  a  result,  the  HURWIN  model  utilizes  data  avail¬ 
able  from  a  specific  storm  to  produce  a  time  history  of  surface  wind  speeds 
and  pressures  over  a  rectangular  grid. 

Within  the  HURWIN  model,  the  total  pressure  which  has  been  depth- 
averaged  over  the  planetary  boundary  layer,  is  prescribed  to  be  a  sum  of  the 
pressure  associated  with  the  tropical  cyclone  which  is  translating  with  the 
storm  at  a  speed  and  a  large-scale  background  pressure  p: 

Ps  =  Pc 


The  well-known  exponential  pressure  law  which  follows  is  used  to  repre¬ 
sent  a  tropical  cyclone  having  a  circularly  symmetric  pressure  field  situated  at 
the  low  pressure  center  of  the  storm: 

Pc  =  Peye  *  ^P^ 


where 

Peye  =  prcssurc  at  the  eye  of  the  storm 
Vp  =  p  -  p  =  the  pressure  anomaly 

R  =  scale  radius,  often  assumed  equivalent  to  the  radius  to  maximum 
wind 

r  =  radial  distance  outwards  from  the  eye  of  the  storm 

The  background  (or  far  field)  pressure  can  be  specified  in  terms  of  a  corre¬ 
sponding  geostrophic  flow  as: 


fkxV=  -(Iv^  (26) 

*  P 

where  fli  is  the  Coriolis  parameter  in  the  vertical  direction  and  Vg  is  the  con¬ 
stant  geostrophic  velocity. 
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The  HURWIN  wind  model  determines  wind  speeds  by  solving  the  equa¬ 
tions  of  horizontal  motion  which  have  been  vertically  averaged  through  the 
depth  of  the  Planetary  Boimdary  Layer  (PBL).  The  two-dimensional  equation 
of  horizontal  motion  is  expressed  in  coordinates  fixed  to  the  Earth: 


-fk  X  V,  =  -(IVp)  (--^)|V, 

at  p  a 


'  \V 

W  I  w 


where 

—  =  total  time  derivative,  given  as  A  +  •  V 

dt  at 

0 

_ =  local  time  change  relative  to  fixed  coordinates 

dt 

V  =  two-dimensional  del  operator 
=  vertically  averaged  horizontal  wind  velocity  vector 
=  mean  air  density 

Kfi  =  horizontal  eddy  viscosity  coefficient  relative  to  the  atmosphere 

CyffQ  =  surface  wind  drag  coefficient 

d  =  depth  of  the  planetary  boundary  layer  (PBL) 

In  the  derivation  of  Equation  27,  the  vertical  advection  of  momentum  is 
assumed  to  be  small  compared  to  horizontal  momentum  advection  and  is  thus 
neglected.  Additionally,  the  shear  stress  at  the  top  of  the  PBL  is  taken  to  be 
zero. 

Substituting  the  pressure  relationships  in  Equations  24  and  26  into  Equa¬ 
tion  27,  the  integrated  wind  field  within  the  PBL  can  be  obtained  from: 


■fk  X  -  n)  -  -  (-VPc)  ^ 


'  -  J  \  ’  yu  -a/  \ - 

dt  ^  p 


d 


The  equations  actually  solved  in  the  HURWIN  wind  model  are  derived  by 
transforming  Equation  28  to  a  moving  Cartesian  coordinate  system  whose 
origin  is  located  at  the  low  pressure  center  of  the  storm. 
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By  neglecting  acceleration  and  horizontal  diffusion  of  momentum  at  the 
edges  of  the  computational  domain,  boundary  conditions  at  these  points  are 
written  as  a  balance  between  the  Coriolis,  pressure  gradient,  and  surface 
frictional  forces. 

The  theoretical  basis  of  the  hurricane  wind  model,  HURWIN,  is  only 
briefly  outlined  here.  A  more  extensive  and  detailed  description  is  given  by 
Cardone,  Greenwood,  and  Greenwood  (1992). 


Numerical  solution  strategy 

The  numerical  discretization  of  the  governing  equations  in  the  HURWIN 
wind  model  is  based  on  a  finite  difference  formulation,  the  details  of  which  are 
given  by  Chow  (1971).  Resolution  of  the  hurricane  is  critical  to  representing 
its  interaction  with  the  environment  (Kurihara  et  al.  1990).  To  best  represent 
the  dynamics  of  the  hurricane  and  maintain  a  computationally  feasible  prob¬ 
lem,  the  discretized  equations  are  solved  over  a  nested  grid  system  consisting 
of  five  rectangular  grids  each  having  uniform  mesh  spacing.  Mesh  spacing 
doubles  for  each  successively  larger  grid  (e.g.,  for  an  innermost  grid  spacing  of 
5  km  (2.7  miles),  remaining  grids  in  the  nested  system  have  mesh  sizes  of 
10  km  (5.4  miles),  20  km  (10.8  miles),  40  km  (21.6  miles),  and  80  km 
(43.2  N  miles)).  The  maximum  areal  extent  of  the  nested  grid  system  corre¬ 
sponds  to  the  areal  coverage  of  the  largest  grid.  The  grid  structure  imple¬ 
mented  within  the  HURWIN  model  is  centered  at  the  eye  of  the  storm  and 
moves  with  the  storm  throughout  the  domain.  This  solution  strategy  provides 
maximum  resolution  near  the  eye  of  the  storm  where  wind  speeds  and  pressure 
gradients  vary  most  rapidly. 

The  HURWIN  model  computations  are  based  on  the  concept  that  a  tropical 
storm  generally  changes  structure  relatively  slowly.  As  such,  the  structure  of  a 
hurricane  travelling  over  open  water  is  well-represented  by  parameters  speci¬ 
fied  at  intervals  of  6  to  24  hr.  These  representative  states  or  “snapshots”  form 
a  portion  of  the  input  to  the  HURWIN  wind  model.  Parameters  specified  at 
6-hr  snapshot  intervals  include:  pressure  at  the  eye  of  the  storm,  far  field  or 
background  pressure,  radius  to  maximum  wind,  stonn  track  direction,  speed  of 
the  storm,  surface  geostrophic  wind  velocity,  direction  of  the  surface  geo- 
strophic  wind,  and  mesh  spacing  of  the  innennost  grid.  The  innermost  nested 
grid  spacing  is  set  at  5  km  for  all  HURWIN  model  simulations  conducted 
herein.  Additional  hurricane  parameters  not  enumerated  above  are  specified 
according  to  average  values  as  suggested  by  Cardone,  Greenwood,  and 
Greenwood  (1992)  and  remain  the  same  for  all  hurricanes  considered. 

Although  the  hurricane  structure  typically  changes  slowly,  storm  position 
can  change  relatively  quickly.  Differences  in  the  time  scale  between  changes 
in  storm  structure  and  position  are  accommodated  by  specifying  the  location  of 
the  storm  center  at  every  hour.  This  hourly  input  to  the  HURWIN  model 
includes  the  latitude  and  longitude  of  the  storm  center  and  the  relative 
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weightings  of  the  hourly  storm  position  between  the  previously  discussed  6-hr 
snapshots  describing  hurricane  parameter  specifications. 

The  6-hr  snapshots  detailing  hurricane  structure  and  the  hourly  data  specify¬ 
ing  storm  position,  both  of  which  are  required  for  operation  of  the  HURWIN 
wind  model,  are  obtained  from  widely  available  meteorological  data  appropri¬ 
ate  to  specific  historical  hurricanes. 

The  parameters  input  into  the  HURWIN  model  at  6-hr  intervals  which 
describe  the  hurricane  wind  field  structure  are  used  to  compute  pressure  gradi¬ 
ent  and  surface  wind  fields  under  steady  state  assumptions.  An  hourly  time 
history  of  the  surface  wind  field  is  then  obtained  through  a  weighted  linear 
interpolation  between  successive  6-hr  computed  wind  fields. 

Despite  the  more  realistic  theoretical  basis  of  the  HURWIN  model,  several 
limitations  arise  in  the  model's  ability  to  represent  hurricane  wind  fields.  The 
most  prominent  shortcoming  relates  to  the  lack  of  dissipation  in  the  wind  field 
of  a  hurricane  making  landfall.  As  a  hurricane  approaches  land  it  generally 
slows  and  winds  dissipate.  Consequently,  winds  on  the  right  side  of  the  hurri¬ 
cane  are  somewhat  reduced  from  full  strength  as  they  encounter  the  land.  As 
the  winds  rotate  over  the  land,  their  speed  is  further  diminished  by  frictional 
forces  resulting  from  the  roughness  of  the  topography.  Consequently,  winds 
on  the  left  side  of  the  hurricane  coming  off  the  land  are  significantly  less  in 
magnitude  than  winds  coming  onto  land  from  the  right  side  of  the  hurricane. 
The  version  of  the  HURWIN  wind  model  used  herein  does  not  represent  the 
processes  associated  with  a  landfalling  hurricane.  Thus,  excessive  winds  near 
the  land  boundary  are  computed  by  the  HURWIN  model  and  may  lead  to  large 
negative  surge  in  coastal  regions  on  the  left  side  of  the  hurricane  and  overpre¬ 
diction  of  the  peak  surge  on  the  right  side  of  the  hurricane. 

Note  that  the  HURWIN  model  does  accommodate  terrains  of  varying 
roughness.  However,  the  incorporation  of  surface  roughness  variations  within 
the  HURWIN  model  does  not  account  for  asymmetry  in  the  winds  of  a  land¬ 
falling  hurricane.  Within  the  HURWIN  model,  the  winds  over  land  can  be 
reduced  by  a  constant  factor  but,  as  the  winds  come  off  the  land  into  the  open 
ocean,  wind  speeds  uniformly  increase  in  accordance  with  the  reduced  fric¬ 
tional  coefficient  specified  over  the  open  ocean.  This  instantaneous,  uniform 
increase  and  subsequent  decrease  of  the  wind  speed  initiated  by  the  specifica¬ 
tion  of  surface  roughness  is  not  physically  realistic.  The  complex  processes 
which  dissipate  winds  in  a  landfalling  hurricane  are  critical  to  generation  of 
storm  surge  at  the  shoreline,  and  yet,  are  not  addressed  by  the  version  of  the 
HURWIN  model  implemented.  For  consistency,  all  surface  roughness  parame¬ 
ters  specified  within  the  HURWIN  wind  model  pertain  to  the  open  ocean. 


Interface  with  ADCIRC-2DDI 

Output  from  the  HURWIN  model  has  been  modified  so  that  computed  wind 
and  pressure  fields  directly  interface  with  the  ADCIRC-2DDI  hydrodynamic 
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model.  The  pressure  field  computed  by  a  HURWIN  model  simulation  is  more 
conveniently  expressed  in  terms  of  an  equivalent  height  of  water  pjpog.  Fur¬ 
thermore,  wind  speeds  generated  by  the  HURWIN  model  are  converted  to 
surface  wind  stresses  using  a  relationship  proposed  by  Garratt  (1977): 


T(i) 


and 


9o 


(29) 


(30) 


where 

\  Tx  =  wind  stresses  in  the  (j)  and  X  directions,  respectively 

PaiJ  Po  =  0.001293  =  ratio  of  air  density  to  average  density  of  seawater 
(Pickard  and  Emery  1982) 

Cp  =  frictional  drag  coefficient  computed  as 
0.001  (0.75  +  0.067 1 Wj) 

|W|  =  magnitude  of  the  wind  velocity 

Wx  =  components  of  the  wind  velocity  vector  in  the  ^  and  X 
directions,  respectively 


The  wind  stress  relationship  of  Garratt  (1977)  has  been  used  successfully 
by  others,  including  Cooper  and  Thompson  (1989);  Westerink,  Stolzenbach, 
and  Connor  (1989);  and  Hearn  and  Holloway  (1990). 


Through  a  bilinear  interpolation,  both  wind  stresses  and  PliPy^g)  values  are 
transferred  directly  from  the  HURWIN  model  nested  rectangular  grids  onto  the 
finite  element  computational  grid  used  by  the  ADCIRC-2DDI  model.  The 
hourly  time  series  of  spatial  wind  stress  and  PJPog  fields  associated  with  a 
specified  hurricane  serve  as  the  meteorological  forcing  for  the  ADCIRC-2DDI 
hydrodynamic  model  in  all  storm  surge  simulations  presented  herein.  Within 
the  ADCIRC-2DDI  model,  linear  interpolation  between  hourly  wind  and  pres¬ 
sure  forcings  yields  meteorological  forcing  appropriate  at  each  time-step 
interval. 
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Hurricanes  of  Record 


Hurricane  Kate 

Hurricane  Kate  is  one  of  the  historical  storms  selected  as  the  meteorological 
forcing  for  the  hydrodynamic  model.  The  track  of  Hurricane  Kate  through  the 
western  North  Atlantic  Ocean  and  the  Gulf  of  Mexico  is  shown  at  6-hr  incre¬ 
ments  in  Figure  1.  Beginning  in  the  deep  western  North  Atlantic  Ocean,  Hur¬ 
ricane  Kate  was  recognized  as  a  tropical  storm  at  18:00  Greenwich  mean  time 
(GMT)  15  November  1985,  and  was  upgraded  24  hr  later  to  a  hurricane.  As 
Kate  moved  northeastward  through  the  Gulf  of  Mexico,  a  minimum  pressure 
of  953  mb  and  a  maximum  sustained  wind  speed  of  approximately  54  m/s  * 
(104.9  knots)  were  recorded  26.5  hr  prior  to  landfall.  At  22:30  GMT 
21  November  1985,  Hurricane  Kate  made  landfall  near  Panama  City,  FL.  The 
tracking  of  Hurricane  Kate  then  continued  through  18:00  GMT  23  November 
1985  until  the  winds  had  significantly  dissipated  and  Kate  was  classified  as  an 
extratropical  storm.  Extensive  meteorological  analysis  of  Hurricane  Kate  as 
well  as  measured  storm  surge  hydrographs  at  10  stations  along  the  Florida 
coast  in  the  vicinity  of  Panama  City  are  reported  by  Garcia  and  Hegge  (1987). 

Simulation  of  Hurricane  Kate  using  the  HURWIN  wind  model  begins  at 
18:00  GMT  15  November  1985  and  ends  8  days  later  at  18:00  GMT 
23  November  1985,  a  duration  of  193  hr.  During  this  period  the  forward 
speed  of  Hurricane  Kate  ranges  from  1.5  m/s'*  (2.9  knots)  to  13.4  m/s'* 

(26.0  knots)  and  the  radius  to  maximum  wind  varies  between  80  km 
(43.2  miles)  and  17  km  (9.2  miles).  A  constant  background  pressure  of 
1,013  mb  is  assumed.  All  of  the  parameters  pertaining  to  Hurricane  Kate  and 
used  within  a  HURWIN  model  simulation  are  obtained  from  the  HURDAT 
tape  (Jarvinen,  Neumaim,  and  Davis  1984;  1993). 


Hurricane  Camille 

Hurricane  CamiUe  is  the  second  historical  storm  used  as  meteorological 
forcing  in  this  study  of  storm  surge  generation.  The  track  of  Hurricane 
CamiUe  through  the  western  Caribbean  Sea  and  the  Gulf  of  Mexico  is  shown 
at  6-hr  increments  in  Figure  2.  As  of  3  August  1969,  Camille  was  tracked  as 
a  tropical  cyclone  having  low  pressures  over  the  Cape  Verde  Islands  in  the 
eastern  Atlantic  Ocean.  During  the  period  covering  12-14  August  1969,  Hur¬ 
ricane  CamiUe  began  its  formation  in  the  western  Caribbean  Sea.  After 
CamiUe  passed  over  Cuba  and  was  well  into  the  Gulf  of  Mexico,  a  minimum 
pressure  of  905  mb  and  a  maximum  sustained  wind  speed  of  about  62  m/s'* 
(120.4  knots)  were  recorded  at  0:00  GMT  17  August  1969.  During  the  same 
period,  water  heights  of  7.5  m  (24.6  ft)  above  mean  sea  level  were  recorded  at 
New  Orleans,  LA.  On  18  August  1969,  Hurricane  CamiUe  made  landfaU  near 
Biloxi,  MS.  The  tracking  of  Hurricane  CamiUe  then  continued  through 
12:00  GMT  22  August  1969  until  the  winds  had  significantly  dissipated. 
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Figure  1 .  Track  of  Hurricane  Kate  through  the  western  North  Atlantic  Ocean 
and  Gulf  of  Mexico  (18:00  GMT  15  November  1985  to  18:00  GMT 
23  November  1985)  shown  at  6-hr  intervals 

Extensive  meteorological  analysis  of  the  formation  of  Hurricane  Camille,  pho¬ 
tos  of  the  storm  track,  recorded  parameters,  and  tidal  heights  were  compiled  by 
Simpson  and  Riehl  (1981). 

Simulations  of  Hurricane  Camille  using  the  HURWIN  wind  model  began  at 
18:(K)  GMT  14  August  1969  and  ended  7-3/4  days  later  at  12:00  GMT 
22  August  1969,  a  duration  of  187  hr.  During  this  period  the  forward  speed  of 
Hurricane  Camille  ranged  from  3.1  m/s'^  to  20.6  m/s'^  (6.0  to  40.0  knots)  and 
the  radius  to  maximum  wind  varied  between  80  km  and  17  km  (43.2  and 
9.2  n.m.).  Throughout  the  simulation,  a  constant  background  pressure  of 
1,013  mb  is  assumed.  All  of  the  parameters  pertaining  to  Hurricane  Camille 
and  used  within  the  HURWIN  model  are  obtained  from  the  HURDAT  tape 
(Jarvinen,  Neumann,  and  Davis  1984;  1993). 
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Figure  2.  Track  of  Hurricane  Camille  through  the  western  Caribbean  Sea 
and  Gulf  of  Mexico  (18:00  GMT  14  August  1969  to  18:00  GMT 
22  August  1969)  shown  at  6-hr  intervals 


Synthetic  Hurricanes 

A  series  of  four  synthetic  hurricanes  are  created  as  meteorological  forcing 
for  the  grid  sensitivity  study  which  investigates  the  relationships  between  grid 
structure,  storm  surge  generation,  and  selected  hurricane  parameters.  The 
synthetic  hurricanes  are  constructed  such  that  characteristics  relating  to  path, 
spatial  scale,  and  forward  velocity  are  isolated  for  study. 

Two  different  hurricane  paths  are  considered,  one  perpendicular  (Path  1) 
and  one  parallel  (Path  3)  to  the  coastline  as  shown  in  Figure  3.  Along  Path  3, 
the  eye  of  the  hurricane  comes  within  60  km  (32.4  n.m.)  of  the  shoreline. 

Two  constant  spatial  scales  of  30  km  (16.2  n.m.)  and  60  km  (32.4  n.m.) 
and  two  constant  storm  velocities  of  15  km/hr  (8.1  knots)  and  25  km/hr 
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(13.5  knots)  are  specified  as 
parameters  for  the  synthetic 
hurricanes.  Each  synthetic 
hurricane  has  a  constant  pres¬ 
sure  deficit  of  80  mb  and  no 
background  surface  winds  are 
applied.  Hurricane  parame¬ 
ters  are  selected  based  on 
averages  of  the  parameters 
associated  with  six  historical 
hurricanes  (Kate,  Eloise, 

Elena,  Hugo,  Danny,  and 
Camille)  and  the  hurricane 
parameters  utilized  in  the 
work  of  Bungpapong,  Reid, 
and  Whitaker  (1985). 

Table  1  summarizes  the 
four  synthetic  hurricanes  and 
their  respective  parameters. 
Hurricane  names  include  three 
numbers  which  identify  the 
parameter  set  associated  with 
each  hurricane.  The  first  number  refers  to  the  specified  pressure  deficit  (0  for 
a  pressure  deficit  of  80  mb).  The  second  number  reflects  the  path  specification 
(1  for  Path  1,  which  is  in  a  direction  perpendicular  to  the  coastline  or  3  for 
Path  3,  which  has  a  segment  parallel  to  the  coastline).  The  third  number  of 
the  hurricane  name  designates  the  combination  of  spatial  scale  and  velocity 
which  applies  to  the  hurricane  (i.e.,  1  refers  to  a  spatial  scale  of  30  km 
(16.2  miles)  and  a  velocity  of  15  km/hr  (8.1  knots)). 


Figure  3. 


Two  synthetic  hurricane  tracks, 
one  perpendicuiar  (Path  1)  and 
one  parallel  to  the  coastline 
(Path  3) 


Table  1 

Characteristics  of  Synthetic  Hurricanes 

Hurricane 

Name 

^m«x 

(km) 

Velocity 

(km/hr) 

Hurricane  Path 
(relative  to 
the  coast) 

Simulation 

Period 

(days) 

Time  of 
Landfall 
(hr) 

H011 

30 

15 

Perpendicular 

8.75 

192 

H012 

60 

15 

Perpendicular 

8.75 

192 

H013 

30 

25 

Perpendicular 

6.00 

126 

H031 

30 

25 

Parallel 

20.50 

— 
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Summary 


Hurricane  wind  and  pressure  fields  are  the  primary  mechanisms  by  which 
stoim  surge  is  generated  in  coastal  regions.  Consequently,  meteorological 
forcing  is  the  driving  force  for  a  storm  surge  hydrodynamic  model.  Realistic 
specification  of  hurricane  wind  and  pressure  forcing  is  a  prerequisite  for  accu¬ 
rate  representation  of  the  storm  surge  response  by  the  hydrodynamic  model. 
Commonly,  empirical  hurricane  models  are  used  but  a  physically  based  conser¬ 
vation  law  model  of  the  planetary  boundary  layer  is  more  appealing  and  per¬ 
haps  more  accurate  in  representing  hurricane  wind  and  pressure  fields.  The 
HURWIN  wind  model  of  Cardone,  Greenwood,  and  Greenwood  (1992)  is  just 
such  a  model  and  is  selected  to  compute  the  meteorological  forcings  for  the 
storm  surge  simulations  herein. 

The  theoretical  basis  and  solution  strategy  of  the  HURWIN  wind  model  are 
presented.  Input  parameters  required  by  the  HURWIN  model  to  generate  wind 
and  pressure  fields  from  specific  historical  hurricanes  are  readily  obtained  from 
available  meteorological  databases,  simplifying  implementation  of  the 
HURWIN  wind  model.  The  interface  between  HURWIN  wind  model  compu¬ 
tations  and  the  wind  and  pressure  forcing  required  by  the  hydrodynamic  model 
ADCIRC-2DDI  is  discussed.  Limitations  in  HURWIN’s  representation  of  a 
landfalling  hurricane  cause  computed  hurricane  wind  fields  which  are 
excessively  high  near  land  and  in  turn  severely  affect  the  accuracy  of  com¬ 
puted  storm  surge  at  the  coastline.  The  effects  of  overpredicted  winds  near 
land  are  in  fact  manifest  in  the  storm  surge  hydrographs  presented  in  subse¬ 
quent  chapters. 

Meteorological  forcing  throughout  this  work  is  obtained  either  from  one  of 
the  historical  hurricanes,  Kate  (November  1985)  or  Camille  (August  1969),  or 
from  one  of  four  synthetic  hurricanes  created  to  isolate  the  influence  of  the 
hurricane  path,  speed,  and  spatial  scale  on  storm  surge  generation.  Details 
regarding  the  characteristics  defining  each  synthetic  storm  are  provided. 
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4  Domain 


Sensitivity 


A  major  shortcoming  in  coastal  ocean  modeling  work  is  the  lack  of  ade¬ 
quate  studies  to  prove  convergence  with  regard  to  grid  structure,  grid  spacing, 
and  domain  size.  In  many  studies  it  is  unclear  whether  the  computed  water 
body  response  is  significantly  aliased  due  to  inadequate  grid  resolution  or  is 
overwhelmed  by  the  imposed  boundary  conditions  and  their  interactions  with 
the  selected  domain.  Only  recently  have  there  been  efforts  to  establish  the 
level  to  which  the  computed  physics  have  converged  by  doing  systematic  grid 
convergence  studies  for  tidal,  wind-driven,  and  large-scale  baroclinic  circula¬ 
tion  (Johns  et  al.  1983a,b;  Le  Provost  and  Vincent  1986;  Bennett  and 
Campbell  1987;  Dietrich,  Roache,  and  Marietta  1990;  Piacsek  and  Allard  1993; 
Dietrich  1993;  Westerink,  Luettich,  and  Muccino  1994;  Luettich  and 
Westerink,  in  preparation).  None  of  these  studies,  however,  specifically 
addresses  the  domain  size  and  its  effect  on  storm  surge  computations. 

The  primary  focus  of  this  chapter  is  to  investigate  the  influence  of  domain 
size  on  hurricane  storm  surge  response.  Domain  size  affects  both  boundary 
condition  specification  and  the  computed  storm  surge  response.  Jarvinen  and 
Lawrence  (1985)  specifically  note  that  the  initial  design  and  location  of  the 
basin  domain  is  a  critical  step  in  the  construction  of  a  storm  surge  model. 
Judicious  selection  of  a  domain  size  results  in  a  computed  storm  surge 
response  that  captures  the  peak  surge  as  well  as  possibly  excited  resonant 
frequencies  documented  by  Bunpapong,  Reid,  and  Whitaker  (1985). 

Storm  surge  modeling  efforts  to  date  have  primarily  used  domains  which 
are  either  limited  to  the  continental  shelf  or  cover  the  shelf  and  portions  of  the 
continental  slope  (e.g.  Flather  1984;  Dendrou,  Moore,  and  Myers  1985; 
Jarvinen  and  Lawrence  1985;  Dube,  Sinha,  and  Roy  1986;  Lardner  and 
Cekirge  1988;  Al-Rabeh,  Eunay,  and  Cekirge  1990;  Hearn  and  Holloway 
1990;  Hubbert,  Leslie,  and  Manton  1990;  Signorini  and  Miller  1992).  Appre¬ 
ciable  development  of  hurricane  storm  surge  occurs  on  the  continental  shelf, 
precluding  accurate  specification  of  the  open  ocean  boundary  condition  in  shelf 
regions.  Nested  domains  are  often  implemented  to  improve  elevation  specifi¬ 
cations  at  the  open  ocean  boundary  of  shelf  domains  (e.g.  Rather  1984; 
Dendrou,  Moore,  and  Myers  1985;  Hubbert,  Leslie,  and  Manton  1990; 
Al-Rabeh,  Eunay,  and  Cekirge  1990).  In  the  nested  domain  approach,  compu¬ 
tations  over  a  large  regional  domain  are  interpolated  onto  the  open  ocean 
boundary  of  a  smaller  shelf  domain.  This  interpolation  from  a  coarsely 
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discretized  domain  to  a  typically  highly  resolved  domain  introduces  numerical 
noise  at  inter-domain  boundaries  (i.e..  Rather  1984;  Hubbert,  Leslie,  and 
Manton  1990)  causing  a  misrepresentation  of  the  dynamics  of  storm  surge 
generation  at  the  open  ocean  boundary  of  the  smaller  domain.  An  extensive 
examination  of  the  storm  surge  response  computed  over  a  shelf  domain  is 
included  as  part  of  the  comprehensive  study  investigating  the  influence  of 
domain  size  on  storm  surge  computations. 

For  this  domain  size  sensitivity  study,  an  actual  hurricane  which  made 
landfall  on  the  Rorida  shelf  in  the  Gulf  of  Mexico  is  applied  over  three 
domains.  Each  domain  is  of  a  different  size  fanning  out  from  the  point  of 
hurricane  landfall  and  covering  successively  larger  regions.  The  first  domain 
is  a  relatively  small  coastal  domain  which  extends  mainly  over  the  continental 
shelf.  The  second  domain  covers  the  entire  Gulf  of  Mexico  basin.  The  final 
domain  extends  well  into  the  deep  Atlantic  Ocean.  Two  different  elevation 
boundary  conditions  are  considered  in  the  simulations,  a  still-water  condition 
and  a  condition  which  incorporates  the  atmospheric  pressure  component  of  the 
meteorological  forcing.  Comparisons  are  made  between  storm  surge  elevations 
computed  over  all  three  domains  using  both  boundary  conditions  to  determine 
the  influence  of  domain  size  and  the  sensitivity  to  boundary  condition  specifi¬ 
cation  on  storm  surge  response.  Furthermore,  computed  storm  surge  elevations 
are  compared  to  measured  water  levels  from  a  historical  hurricane.  Results  of 
this  domain  size  sensitivity  study  provide  criteria  for  domain  size  selection 
with  regard  to  the  accurate  generation  and  propagation  of  hurricane  storm 
surge. 


Hydrodynamic  Domain  Descriptions 

Three  hydrodynamic  domains  of  widely  varied  sizes  are  considered  in  this 
investigation.  The  domain  sizes  are  selected  to  clearly  demonstrate  the  rela¬ 
tionship  between  the  domain  and  open  boundary  elevation  specification.  The 
two  smaller  domains  are  constructed  to  correspond  with  domains  used  in 
recent  storm  surge  modeling  efforts.  Hurricane  Kate  is  the  meteorological 
forcing  used  for  the  domain  size  comparison  simulations.  Thus,  all  domains 
used  here  expand  outward  from  the  landfall  region  of  Hurricane  Kate  near 
Panama  City,  FL,  and  cover  increasingly  larger  areas.  Relative  sizes  of  the 
three  domains  considered  are  shown  in  Figure  4  and  the  exact  areal  extent  of 
each  domain,  as  well  as  other  characteristics  pertinent  to  the  domains  and  their 
discretizations,  are  summarized  in  Table  2. 

The  smallest  domain  considered  is  the  Rorida  coast  domain  pictured  in 
Figure  4a.  The  Rorida  coast  domain  is  a  semicircular  basin  similar  to  the  one 
used  by  the  National  Weather  Service  at  Pensacola  Bay  in  conjunction  with  the 
SLOSH  storm  surge  model  (Jelesnianski  1979;  Jarvinen  and  Lawrence  1985; 
Shaffer,  Jelesnianski,  and  Chen  1986;  Jelesnianski,  Chen,  and  Shaffer  1992; 
Chen,  Shaffer,  and  Kim  1993).  The  Rorida  coast  domain  extends  radially 
outward  into  the  Gulf  of  Mexico  from  the  shoreline  surrounding  Panama  City, 
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Figure  4.  Boundaries  of  the  (a)  Florida  Coast  domain,  (b)  Gulf  of  Mexico 
domain,  and  (c)  east  coast  domain 


Table  2 

Characteristics  of  Model  Domains 

Domain  Name 

Area 

(km*) 

Max 

Depth 

(m) 

Discretization 

Grid  Size  (km) 

Nodes 

Elements 

Max 

Min 

Florida  coast 

5.07  X  10“ 

1,094 

1,451 

2,326 

32,5 

0.5 

Gulf  of  Mexico 

1,41  X  10® 

3,781 

6,325 

11,441 

50,0 

0,5 

East  coast 

8.35  X  10® 

7,765 

22,71 1 

41,709 

98.3 

0.5 

FL,  with  a  radius  of  approximately  175  km  (109  miles).  The  outer  arc  of  the 
semicircular  basin  forms  the  open  ocean  boundary  of  the  domain.  Bathymetric 
data  for  the  Florida  coast  domain  shown  in  Figure  5  are  obtained  from  the 
topographic  database  ET0P05  from  the  National  Center  for  Atmospheric 
Research  and  supplemented  by  the  National  Oceanic  and  Atmospheric  Admini¬ 
stration  (NOAA)  Digital  U.S.  Coastal  Hydrography  sounding  database  (distrib¬ 
uted  by  NOAA  National  Geophysical  and  Solar-Terrestrial  Data  Center  in 
Boulder,  CO)  along  the  Florida  shelf  and  coast.  A  minimum  depth  of  3  m 
(10  ft)  is  imposed  throughout  the  Florida  coast  domain  to  eliminate  complete 
drying  in  computational  elements  along  the  shoreline  and  in  embayments.  A 
maximum  depth  of  nearly  1,100  m  (3,6(X)  ft)  is  recorded  at  the  outer  limit  of 
the  Florida  coast  domain,  but  as  is  seen  in  Figure  5,  much  of  the  Florida  coast 
domain  lies  on  the  continental  shelf  at  depths  less  than  130  m  (426  ft). 
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Discretization  of  the 
Florida  coast  domain,  shown 
in  Figure  6,  is  achieved  using 
1,451  nodes  and  2,326  ele¬ 
ments.  The  computational 
grid  for  the  Florida  coast 
domain  is  extremely  detailed 
along  the  shoreline  and  has 
considerable  refinement  in  the 
near  coastal  region  resulting 
in  a  minimum  node-to-node 
spacing  of  approximately 
0.5  km  (0.3  mile).  Maximum 
node-to-node  spacing  of 
32.5  km  (20.2  statute  miles) 
is  found  in  the  deeper  regions 
of  the  Florida  coast  domain. 


BaafiSiBiSgS&esfeai 


Figure  5.  Bathymetry  contours  for  the  Florida  xhe  second  domain  con- 
coast  domain  sidered,  the  Gulf  of  Mexico 

domain  shown  in  Figure  4b, 
includes  the  entire  Gulf  of 
Mexico  and  is  similar  to 
domains  used  by  several  other 
investigators  for  storm  surge 
modeling  studies  within  the 
Gulf  of  Mexico  (e.g. 
Bunpapong,  Reid,  and 
Whitaker  1985;  Westerink  et 
al.  1992a).  The  Gulf  of 
Mexico  domain  is  comprised 
of  the  Florida  coast  domain 
and  all  surrounding  regions  in 
the  Gulf  of  Mexico.  One 
open  ocean  boundary  is 
located  across  the  Strait  of 
Florida  (approximately  from 
Cape  Sable,  Florida,  to 

Figure  6.  Florida  coast  domain  discretization  Havana,  Cuba)  and  another  is 
located  across  the  Yucatan  Channel  (approximately  from  Cancun,  Mexico,  to 
Cabo  San  Antonio,  Cuba).  Bathymetry  in  the  Gulf  of  Mexico  domain, 
depicted  in  Figure  7,  is  again  taken  from  the  ET0P05  database  and  in  regions 
along  the  Florida  coast  and  shelf  it  is  supplemented  by  the  NOAA  Digital 
U.S.  Coastal  Hydrography  sounding  database.  In  the  portions  of  the  Gulf  of 
Mexico  domain  which  correspond  to  the  Florida  coast  domain,  the  topography 
is  identical  to  that  specified  for  the  Florida  coast  domain.  The  range  of  depths 
in  the  Gulf  of  Mexico  domain  vary  between  a  minimum  of  3  m  (10  fit)  and  a 
maximum  just  imder  3,800  m  (12,500  ft). 
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The  discretization  of  the 
Gulf  of  Mexico  domain  is 
presented  in  Figure  8  and  is 
comprised  of  6,325  nodes  and 
11,441  elements.  In  the 
region  encompassing  the 
Rorida  coast  domain,  the 
discretization  of  the  Gulf  of 
Mexico  domain  is  the  same  as 
that  in  the  Florida  coast 
domain.  Nodal  spacings 
within  the  Gulf  of  Mexico 
domain  discretization  range 
from  a  minimum  of  0.5  km 
(0.3  statute  miles)  along  the 
shoreline  and  in  the  coastal 
region  of  western  Florida  to  a 
maximum  of  50  km  (31.1  sta¬ 
tute  miles)  in  the  deepest 
waters  of  the  Gulf  of  Mexico. 


The  final  and  largest 
domain  studied,  shown  in 
Figure  4c,  is  the  east  coast 
domain,  which  has  been  pre¬ 
viously  used  by  Westerink, 
Luettich,  and  Muccino  (1994) 
to  study  tides  in  the  western 
North  Atlantic.  The  east 
coast  domain  encompasses  the 
western  North  Atlantic  Ocean, 
the  Caribbean  Sea,  and  the 
Gulf  of  Mexico  and  has  been 
constmcted  such  that  both  the 
Gulf  of  Mexico  and  Florida 
coast  domains  are  contained 
within  the  east  coast  domain. 
A  single  deep  Atlantic  Ocean 
boundary  within  the  east  coast 
domain  extends  from  Glace  Bay,  Nova  Scotia,  to  the  vicinity  of  Corocora 
Island  in  eastern  Venezuela  along  the  60  °W  meridian.  All  other  boundaries 
are  defined  by  the  eastern  coastlines  of  North,  Central,  and  South  America. 
Topography  within  the  domain,  depicted  in  Figure  9,  includes  the  continental 
shelf,  whose  depths  range  from  an  imposed  minimum  between  3  m  and  7  m 
(10  ft  and  426  ft)  to  130  m  (426  ft)  at  the  shelf  break,  the  continental  slope 
which  has  a  typical  depth  range  of  130  m  to  3,000  m  (426  ft  to  9,800  ft),  and 
the  continental  rise  and  deep  ocean  where  depths  increase  upwards  from 
3,000  m  to  almost  8,000  m  (9,800  ft  to  26,000  ft).  Depths  within  the  Gulf  of 
Mexico  region  are  the  same  as  those  specified  for  the  Gulf  of  Mexico  domain. 


Figure  7.  Bathymetry  contours  for  the  Gulf  of 
Mexico  domain 
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Figure  8.  Gulf  of  Mexico  domain  discretization 


Chapter  4  Domain  Size  Sensitivity 


Bathymetry  in  regions  outside 
the  Gulf  of  Mexico  is 
obtained  from  the  ETOP05 
database. 

The  discretization  of  the 
east  coast  domain,  shown  in 
Figure  10,  contains 
22,711  nodes  and  41,709  ele¬ 
ments,  reasonable  numbers 
considering  the  level  of 
refinement  and  the  areal 
extent  of  the  domain.  The 
discretization  of  the  Gulf  of 
Mexico  domain  is  exactly 
incorporated  into  the  east 
coast  domain  discretization. 
The  inclusion  of  deep  Atlantic 
ocean  regions  within  the  east 
coast  domain  extends  the 
range  of  grid  spacings  con¬ 
tained  within  the  east  coast 
domain  discretization.  Mini¬ 
mum  node-to-node  spacing  is 
approximately  0.5  km 
(0.3  statute  miles)  along  the 
Florida  shoreline,  while  in  the 
deep  Atlantic  Ocean,  spacing 
increases  to  about  98  km 
(60. 1  miles)  or  approximately 
1.4  °.  A  large  computational 
domain  such  as  the  east  coast 
domain  is  manageable 
because  of  an  optimal  grid- 
ding  strategy  and  the  flexi¬ 
bility  inherent  in  the  finite 
element  method. 


Comparison  of  Storm  Surge  Predictions  Using 
Three  Domain  Sizes 

Simulations  were  conducted  to  investigate  the  influence  of  domain  size  and 
boundary  condition  specification  on  the  storm  surge  and  resonant  mode  genera¬ 
tion  associated  with  Hurricane  Kate.  Domain  size  is  considered  using  the  three 
previously  described  domains.  Boundary  condition  specification  is  examined 


Figure  9.  Bathymetry  contours  for  the  east 
coast  domain 


Figure  10.  East  coast  domain  discretization 
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by  comparing  a  still-water  boundary  condition  where  water  elevation  is  set 
equal  to  mean  sea  level  with  a  boundary  condition  that  partially  accounts  for 
meteorological  forcing  at  the  open  ocean  boundary  by  imposing  an  inverted 
barometer  effect.  An  inverted  barometer  pjp^g  is  simply  the  height  to  which 
seawater  will  rise  due  to  pure  static  pressure  forcing.  For  this  series  of  simula¬ 
tions,  wind  stress  and  pressure  forcing  from  Hurricane  Kate  are  applied  at  the 
boundary  and  on  the  interior  of  the  domain,  while  tidal  forcing  is  neglected 
both  on  the  interior  and  at  the  open  ocean  boundaries. 

Model  parameters  are  identical  for  all  simulations  so  that  results  for  differ¬ 
ent  domains  using  each  boundary  condition  can  be  compared.  The  convective 
and  finite  amplitude  terms  are  not  included  in  the  governing  equations  for 
these  simulations  due  to  instabilities  caused  by  near  drying  elements.  Conse¬ 
quently,  the  only  nonlinear  term  included  in  the  governing  equations  is  bottom 
friction.  The  bottom  friction  coefficient  is  constant  and  equ^  to  0.003  over  all 
domains.  The  GWCE  parameter  Tq,  which  represents  the  balance  between  the 
primitive  continuity  and  wave  equation  portion  of  the  GWCE,  is  defined  equal 
to  0.001  (Kolar  et  al.  1994b). 

Simulations  are  spun  up  from  homogeneous  initial  conditions  using  a  1-day 
ramp  in  time.  Application  of  the  hyperbolic  ramp  function  reduces  the  excita¬ 
tion  of  nonphysical  short  wavelength  frequencies.  An  identical  ramp  function 
of  1  day  is  applied  to  the  wind  and  pressure  forcing  as  well  as  the  inverted 
barometer  boundary  condition  when  used.  Actual  simulations  began  at 
12:00  GMT  15  November  1985  and  ran  over  8.25  days  (including  a  1-day 
ramp-up  period).  During  the  first  6  hr  of  the  simulation,  the  initial  hurricane 
wind  and  pressure  forcings  were  held  stationary.  Thereafter,  storm  surge  com¬ 
putations  used  the  time-varying  wind  and  pressure  fields.  A  time-step  of 
45  sec  is  used  throughout  the  simulation  period.  No  calibration  or  tuning  of 
parameters  is  performed  in  either  the  weather  model  or  in  the  hydrodynamic 
model. 

As  part  of  this  domain  size  study,  118  elevation  stations  were  placed 
throughout  the  three  domains.  Ten  elevation  stations  were  concentrated  on  the 
Florida  coast  extending  from  St.  Marks,  Florida,  to  the  Mississippi  Sound. 
Other  stations  were  placed  along  the  Gulf  of  Mexico  coastline,  at  the  open 
ocean  boundaries  of  each  domain,  and  on  tracks  parallel  to  the  path  of  Hurri¬ 
cane  Kate.  Samplings  of  these  elevation  stations  are  shown  in  Figures  1 1 
and  12. 

Storm  surge  response  at  nine  of  the  stations  along  the  Florida  coast  (St. 
Marks,  Shell  Point,  Turkey  Point,  Carrabelle,  Apalachicola,  Panama  City, 
Alligator  Bayou,  Destin,  and  Pensacola)  and  Station  T6.3  on  the  continental 
shelf,  is  shown  in  Figures  13-22.  The  computed  storm  surge  hydrographs  in 
Figures  13a-22a  were  obtained  using  a  still- water  boundary  condition  and  in 
Figures  13b-22b  using  an  inverted  barometer  boundary  condition.  The  storm 
surge  profiles  at  St.  Marks,  Shell  Point,  Turicey  Point,  Carrabelle,  and 
Apalachicola  (Figures  13-17)  are  representative  of  conditions  on  the  right-hand 
side  of  the  hurricane.  Storm  surge  elevations  at  Alligator  Bayou,  Destin,  and 
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Figure  13.  Computed  storm  surge  from  Hurricane  Kate  using  (a)  a  still-water  boundary  condi¬ 
tion,  and  (b)  an  inverted  barometer  boundary  condition  at  St.  Marks,  FL 
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Figure  14.  Computed  storm  surge  from  Hurricane  Kate  using  (a)  a  still-water  boundary  condi¬ 
tion,  and  (b)  an  inverted  barometer  boundary  condition  at  Shell  Point,  FL 
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Figure  15.  Computed  storm  surge  from  Hurricane  Kate  using  (a)  a  still-water  boundary  condi¬ 
tion,  and  (b)  an  inverted  barometer  boundary  condition  at  Turkey  Point.  FL 
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Figure  16.  Computed  storm  surge  from  Hurricane  Kate  using  (a)  a  still-water  boundary  condi¬ 
tion,  and  (b)  an  inverted  barometer  boundary  condition  at  Carrabelle,  FL 
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Figure  17.  Computed  storm  surge  from  Hurricane  Kate  using  (a)  a  still-water  boundary  condi¬ 
tion,  and  (b)  an  inverted  barometer  boundary  condition  at  Apalachicola,  FL 
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Figure  18.  Computed  storm  surge  from  Hurricane  Kate  using  (a)  a  still-water  boundary  condi¬ 
tion,  and  (b)  an  inverted  barometer  boundary  condition  at  Panama  City,  FL 
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Figure  19.  Computed  storm  surge  from  Hurricane  Kate  using  (a)  a  still-water  boundary  condi¬ 
tion,  and  (b)  an  inverted  barometer  boundary  condition  at  Alligator  Bayou,  FL 
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Figure  20.  Computed  storm  surge  from  Hurricane  Kate  using  (a)  a  still-water  boundary  condi¬ 
tion,  and  (b)  an  inverted  barometer  boundary  condition  at  Destin,  FL 
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3.5 


Hurricane  Kate 


Figure  21 .  Computed  storm  surge  from  Hurricane  Kate  using  (a)  a  still-water  boundary  condi¬ 
tion,  and  (b)  an  inverted  barometer  boundary  condition  at  Pensacola,  FL 
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Figure  22.  Computed  storm  surge  from  Hurricane  Kate  using  (a)  a  still-water  boundary  condi¬ 
tion,  and  (b)  an  inverted  barometer  boundary  condition  at  Station  T6.3 
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Pensacola  (Figures  19-21)  exhibit  behavior  characteristic  of  the  left-hand  side 
of  the  hurricane.  The  Panama  City  storm  surge  hydrograph  (Figure  18)  depicts 
water  elevations  at  the  point  of  landfall  of  Hurricane  Kate.  A  representative 
storm  surge  profile  on  the  shelf  is  shown  in  Figure  22  for  Station  T6.3  located 
on  the  continental  shelf  approximately  15  km  from  the  shelf  break  at  a  depth 
of  nearly  130  m. 

In  aU  of  the  storm  surge  hydrographs  in  Figures  13-22,  the  peak  surge  for 
the  Gulf  of  Mexico  and  east  coast  domains  closely  correspond,  whereas  the 
elevation  response  computed  over  the  Florida  coast  domain  is  significantly 
smaller,  particularly  on  the  right-hand  side  of  the  hurricane.  For  a  given  boun¬ 
dary  condition  specification,  similar  oscillatory  patterns  or  modes  are  excited 
in  the  storm  surge  response  at  all  10  stations  shown  in  Figures  13-22.  How¬ 
ever,  in  comparing  the  computed  storm  surge  response  between  the  east  coast 
and  Gulf  of  Mexico  domains,  different  oscillatory  behavior  is  evident.  For  the 
Gulf  of  Mexico  domain,  resonant  modes  are  more  pronounced  when  a  still- 
water  boundary  condition  is  used  (Figures  13a-22a)  than  when  the  inverted 
barometer  boundary  condition  is  implemented  (Figures  13b-22b).  The  modes 
in  the  east  coast  domain  solution  in  Figures  13-22  exhibit  no  sensitivity  to 
boxmdary  condition  specification.  Since  aU  domains  have  identical  discretiza¬ 
tions  over  corresponding  regions,  and  simulations  were  conducted  using  identi¬ 
cal  model  parameters  and  wind  and  pressure  forcing,  differences  between  the 
model  responses  are  due  solely  to  the  domain  size  and/or  the  boundary  condi¬ 
tion  specification. 

To  investigate  the  origin  of  differences  in  the  storm  surge  response  over  the 
three  domains,  elevation  contours  showing  the  progression  of  Hurricane  Kate 
through  the  Gulf  of  Mexico  toward  Panama  City,  Florida,  are  shown  in  Fig¬ 
ure  23  for  successive  time  periods.  At  20;(X)  GMT  20  November  1985  (Fig¬ 
ure  23a),  a  widespread  inverted  barometer  centered  at  the  hurricane  eye  is 
evident.  Already  water  is  accumulating  and  elevation  is  increasing  in  the 
northeastern  Gulf  of  Mexico,  especially  on  the  northeast  Florida  shelf.  Ten 
hours  later  at  6:00  GMT  21  November  1985  (Figure  23b),  as  the  eye  of  the 
hurricane  approaches  the  Florida  shelf,  the  inverted  barometer  effect  has  dimin¬ 
ished  and  its  structure  is  deteriorating.  High  water  elevations  on  the  shelf, 
particularly  in  shallow  coastal  areas,  indicate  that  hurricane  winds  have  pushed 
water  up  on  the  continental  shelf.  The  last  snapshot,  shown  in  Figure  23c,  is 
taken  8  hr  prior  to  landfall  of  the  hurricane  at  18:00  GMT  21  November  1985. 
Hurricane  Kate  is  now  on  the  continental  shelf  and  the  surface  elevation 
response  no  longer  exhibits  any  inverted  barometer  structure.  Water  continues 
to  pile  up  on  the  shelf  in  a  pronounced  way  on  the  right-hand  side  of  the 
hurricane. 

The  Florida  coast  domain  is  located  almost  entirely  within  the  region  of 
intense  storm  surge  generation  shown  in  Figure  23.  In  particular,  the  cross¬ 
shelf  boundaries  of  the  Florida  coast  domain  are  located  in  an  area  of  signifi¬ 
cant  surge  where  hurricane  winds  have  driven  water  up  onto  the  shelf.  Even 
an  inverted  barometer  forcing  specified  at  these  cross-shelf  boundaries  signifi¬ 
cantly  underestimates  the  storm  surge  elevations  that  physically  occur  at  these 
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Figure  23.  Storm  surge  elevations  (in  cm)  in  the  northeast  Gulf  of  Mexico  due 
to  Hurricane  Kate  on  (a)  20:00  GMT  20  November  1985, 

(b)  6:00  GMT  21  November  1985,  and  (c)  18:00  GMT 
21  November  1985 
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boundaries  and  are  represented  using  the  Gulf  of  Mexico  and  east  coast 
domains.  The  storm  surge  hydrograph  at  Station  T8.5  (Figure  24),  which  is 
located  at  the  cross-shelf  boundary  of  the  Florida  coast  domain  at  a  depth  of 
approximately  80  m  (262  ft),  clearly  demonstrates  the  inadequacy  of  the 
inverted  barometer  specified  at  these  boundaries.  The  inverted  barometer 
elevations  are  significantly  less  than  storm  surge  elevations  computed  over  the 
east  coast  and  Gulf  of  Mexico  domains  at  the  same  point.  Since  the  storm 
surge  generated  on  the  shelf  by  hurricane  winds  is  not  known  a  priori,  an 
appropriate  elevation  boundary  condition  cannot  be  specified  at  cross-shelf 
boundaries.  Furthermore,  these  cross-shelf  boundaries  comprise  nearly  the 
entire  open-ocean  boundary  of  the  small  shelf  domain.  Consequently,  a 
domain  such  as  the  Florida  coast  domain,  which  is  small  in  size  relative  to  the 
spatial  scale  of  a  hurricane  and  is  located  primarily  on  the  continental  shelf, 
cannot  be  used  to  obtain  an  accurate  and  physically  relevant  storm  surge 
response. 


Figure  24.  Computed  storm  surge  from  Hurricane  Kate  using  an  inverted 
barometer  boundary  condition  at  Station  T8.5 


The  oscillatory  behavior  of  the  storm  surge  elevations  or  surge  forerunner 
effect  computed  over  the  east  coast  and  Gulf  of  Mexico  domains  exhibited  in 
Figures  13-22  can  be  attributed  to  basin-wide  resonant  modes  which  exist  in 
the  Gulf  of  Mexico  as  documented  by  Platzman  (1972);  Reid  and  Whitaker 
(1981);  and  Bunpapong,  Reid,  and  Whitaker  (1985).  These  modes  are  not 
excited  in  the  Florida  coast  domain  due  to  the  small  size  of  the  domain  and  the 
omission  of  the  shelf  from  adjacent  basin  interaction  in  the  specified  boundary 
forcing.  On  the  contrary,  resonant  modes  are  excited  in  the  Gulf  of  Mexico 
domain.  In  fact,  open  boundary  elevations  specified  at  entrances  to  the  Gulf  of 
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Mexico  significantly  influence  the  setup  of  resonant  modes  in  the  Gulf.  This 
is  demonstrated  in  Figures  13-22  for  the  Gulf  of  Mexico  domain,  where  oscil¬ 
lations  in  storm  surge  elevation  are  more  pronounced  when  using  a  still-water 
boundary  condition  (Figures  13a-22a)  as  opposed  to  an  inverted  barometer 
boundary  condition  (Figures  13b-22b).  Platzman  (1972)  and  Reid  and 
Whitaker  (1981)  also  found  that  the  frequency  of  modes  in  the  Gulf  of  Mexico 
varied  with  the  application  of  different  boundary  conditions.  One  concludes, 
then,  that  the  resonant  modes  in  the  Gulf  of  Mexico  may  be  difficult  to  model 
with  the  Gulf  of  Mexico  domain  due  to  the  sensitivity  of  these  modes  to 
boundary  condition  specification  and  the  associated  shortcomings  in  the  repre¬ 
sentation  of  basin-to-basin  dynamics. 

A  comparison  of  water  level  response  at  all  coastal  stations  using  two  dif¬ 
ferent  elevation  conditions  along  the  east  coast  domain  open-ocean  boundary 
(e.g.  Figure  25  at  Apalachicola,  Florida)  clearly  demonstrates  that  the  influence 
of  boundary  condition  specification  is  minimized  when  using  the  east  coast 
domain.  This  insensitivity  to  the  specification  of  boundary  conditions  in  the 
deep  Atlantic  as  well  as  the  fact  that  basin-to-basin  interactions  are  more  faith¬ 
fully  represented  suggests  that  the  flow  physics  at  this  boundary  need  not  be 
precisely  known  and  resonant  modes  in  the  Gulf  of  Mexico  are  properly  set  up 
in  the  east  coast  domain. 


Figure  25.  Comparison  of  storm  surge  computed  over  the  east  coast  domain 
using  both  still-water  boundary  condition  and  an  inverted  barome¬ 
ter  boundary  condition  at  Apalachicola.  FL 


Thus,  when  using  the  east  coast  domain,  hurricanes  progress  through  the 
domain  in  a  realistic  fashion  causing  water  to  properly  accumulate  on  the 
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continental  shelf  for  the  following  three  reasons:  (a)  factors  affecting  storm 
surge  generation  are  more  fuUy  incorporated  (i.e.,  bathymetry,  shelf  to  adjacent 
basin  interaction,  basin-  to-basin  interaction),  (b)  resonant  modes  are  more 
realistically  excited  without  relying  on  the  incorporation  of  this  information 
into  the  boundary  conditions,  and  (c)  the  influence  of  the  boundary  condition 
specification  is  generally  minimized.  This  initial  series  of  simulations  clearly 
indicates  that  the  east  coa^t  domain  is  optimal  for  both  the  representation  of 
accurate  primary  storm  surge  and  surge  forerunner  effects. 


Comparison  Between  Observed  and  Computed 
Storm  Surge  Elevations  Over  the  East  Coast 
Domain 

Uncalibrated,  computed  stonn  surge  profiles  at  several  elevation  stations  are 
now  compared  to  measured  storm  surge  data  for  verification  of  the  model 
results.  As  a  consequence  of  the  analysis  in  the  previous  section,  the  east 
coast  domain  is  utilized  for  this  storm  surge  simulation.  As  in  prior  simula¬ 
tions,  wind  stress  and  pressure  forcing  from  Hurricane  Kate  are  applied 
throughout  the  domain.  In  addition,  tidal  potential  forcing  using  five  tidal 
constituents,  K^,  O^,  P^,  M2,  and  Sj,  is  specified  on  the  interior  of  the  domain. 
Table  3  details  the  frequencies,  nodal  factors,  equilibrium  arguments,  and  Earth 
elasticity  factors  used  for  each  tidal  constituent.  The  values  of  these  tidal 
parameters  are  defined  with  reference  to  0:00  GMT  19  November  1985. 

Along  the  open-ocean  boundary,  an  inverted  barometer  boundary  condition  is 
applied  in  combination  with  tidal  elevation  forcing.  The  same  five  constitu¬ 
ents,  K^,  Oj,  Pp  M2,  and  S2,  are  forced  along  the  boundary  with  elevations 
obtained  from  the  results  of  Schwiderski's  global  model  (Schwiderski  1979, 
1981a-d).  However,  the  ampUtude  of  the  M2  tide  on  the  open-ocean  boundary, 
as  computed  by  Schwiderski,  is  reduced  by  10  percent,  a  correction  in  accor¬ 
dance  with  the  error  analysis  of  tidal  model  computations  over  the  western 
North  Atlantic  Ocean  by  Westerink,  Luettich,  and  Muccino  (1994a). 


Table  3 

Tidal  Potential  Constants  Used  for  Storm  Surge  Simulations 

Constituent 

Period 

(hrs) 

Potential 

Amplitude 

(m) 

Earth  Tide 
Reduction 
Factor 

Nodal 

Factor 

Equilibrium 

Argument 

(degrees) 

K,  luni-solar 

23.934470 

0.141565 

0.736 

1.086 

321.3 

0^  principal  lunar 

25.819342 

0.100514 

0.695 

1.140 

271.35 

P,  principal  solar 

24.065890 

0.046843 

0.706 

1,000 

33.00 

M2  principal  lunar 

12.420601 

0.242334 

0.693 

0.974 

230.34 

S2  principal  solar 

12.000000 

0.112841 

0,693 

1.000 

0,00 
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Simulations  used  for  verification  of  the  computed  storm  surge  predictions 
are  conducted  using  a  fully  nonlinear  model  formulation  with  the  exception  of 
the  finite  amplitude  terms.  These  terms  are  not  included  in  the  simulations 
because  of  instabilities  caused  by  near  drying  elements  in  coastal  areas.  All 
other  model  parameters  such  as  the  bottom  friction  coefficient,  the  GWCE 
parameter,  and  the  time-step  remain  as  specified  for  the  domain  comparison 
simulations  detailed  in  the  previous  section.  Again,  all  parameters  in  both  the 
hydrodynamic  and  wind  models  remain  uncalibrated. 

Verification  of  the  tidal  responses  throughout  the  east  coast  domain  is 
undertaken  prior  to  simulation  of  the  storm  surge.  Tides  are  computed  for 
45  days  with  a  ramp-up  period  of  12  days.  Tidal  elevations  are  recorded  after 
30  days  of  computation  to  ensure  dissipation  of  free  modes  within  the  Gulf  of 
Mexico  and  Caribbean  Sea  as  suggested  by  Westerink,  Luettich,  and  Muccino 
(1994).  At  77  stations  placed  throughout  the  east  coast  domain,  computed 
tidal  elevations  are  compared  to  measured  water  levels  obtained  at  these  same 
stations.  Details  regarding  the  location  of  all  77  stations  and  the  source  of 
measured  tidal  data  are  found  in  Westerink,  Luettich,  and  Muccino  (1994). 
Results  are  presented  for  a  sampling  of  six  tidal  stations  whose  locations  are 
given  in  Figures  11  and  12  (i.e..  Key  West,  Cedar  Key,  St.  Maries,  Alligator 
Bayou,  Southwest  Pass,  and  the  outer  Florida  shelf  station).  Tidal  predictions 
shown  in  Figures  26-31  agree  well  with  the  measured  data.  The  slight  down¬ 
ward  adjustment  of  the  M2  amplitude  specified  on  the  open-ocean  boundary 
produces  significantly  improved  tidal  elevations  at  all  stations  as  compared  to 
those  computed  by  Westerink,  Luettich,  and  Muccino  (1994). 

Storm  surge  simulations  are  spun  up  from  homogeneous  initial  conditions 
using  a  1-day  hyperbolic  ramp  in  time.  Actual  simulations  begin  at  0:00  GMT 
10  October  1985  and  run  a  total  of  45  days.  During  the  first  36.75  days  of 
simulation  (including  the  1-day  ramp-up  period),  tidal  forcing  alone  is  applied 
before  the  hurricane  enters  the  domain.  After  this  period  of  pure  tidal  forcing, 
initial  hurricane  wind  and  pressure  forcings  are  held  stationary  for  a  6-hr 
ramp-up  period.  Thereafter,  storm  surge  computations  use  the  time-varying 
wind  and  pressure  forcings  and  proceed  for  8  days  ending  at  18:00  GMT 
23  November  1985.  Storm  surge  elevations  are  recorded  every  3.75  min  at 
eight  elevation  stations  (shown  in  Figure  11)  located  along  the  Florida  coast¬ 
line  where  measured  storm  surge  heights  are  also  available. 

Storm  surge  hydrographs,  shown  in  Figures  32-39,  compare  computed 
numerical  storm  surge  elevations  with  the  measured  water  levels  during  Hurri¬ 
cane  Kate.  Conditions  on  the  right-hand  side  of  the  hurricane,  given  by  the 
computed  storm  surge  profiles  in  Figures  32-36,  correspond  reasonably  well  to 
measured  storm  surge  elevations.  Deviations  from  the  observed  water  levels 
include  a  premature  and  overpredicted  peak  surge  and  slight  underprediction 
following  the  storm  passage.  At  stations  located  at  or  left  of  the  hurricane 
landfall  position  (Figures  37-39),  computed  storm  surge  hydrographs  exhibit  an 
enhanced  and  slightly  elevated  tidal  cycle  while  actual  elevations  show  a  more 
pronounced  surge  which  extends  over  about  1.5  days. 
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Figure  26.  Comparison  of  computed  and  measured  tidal  elevations  over  the 
east  coast  domain  at  Key  West,  FL 


Figure  27.  Comparison  of  computed  and  measured  tidal  elevations  over  the 
east  coast  domain  at  Cedar  Key,  FL 
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NOVEMBER,  1985 


Figure  28.  Comparison  of  computed  and  measured  tidal  elevations  over  the 
east  coast  domain  at  St.  Marks,  FL 


NOVEMBER,  1985 


Figure  29.  Comparison  of  computed  and  measured  tidal  elevations  over  the 
east  coast  domain  at  Alligator  Bayou.  FL 
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Figure  30.  Comparison  of  computed  and  measured  tidal  elevations  over  the 
east  coast  domain  at  Southwest  Pass,  LA 


Figure  31 .  Comparison  of  computed  and  measured  tidal  elevations  over  the 
east  coast  domain  at  the  outer  Florida  shelf  station 
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Hurricane  Kate 


Figure  32.  Comparison  of  computed  and  measured  storm  surge  elevations  for 
Hurricane  Kate  over  the  east  coast  domain  at  St.  Marks,  FL 


NOVEMBER,  1985 


Figure  33.  Comparison  of  computed  and  measured  storm  surge  elevations  for 
Hurricane  Kate  over  the  east  coast  domain  at  Shell  Point,  FL 
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Figure  34.  Comparison  of  computed  and  measured  storm  surge  elevations  for 
Hurricane  Kate  over  the  east  coast  domain  at  Turkey  Point,  FL 


Figure  35.  Comparison  of  computed  and  measured  storm  surge  elevations  for 
Hurricane  Kate  over  the  east  coast  domain  at  Carrabelle,  FL 
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Deviations  from  the  observed  storm  surge,  particularly  the  overprediction  at 
stations  east  of  Panama  City,  Florida,  and  the  lack  of  a  sustained  surge  at 
stations  west  of  Panama  City,  Florida,  may  be  indicative  of  an  artificially  high 
wind  field  forcing  generated  by  the  HURWIN  wind  model.  A  limitation  of  the 
HURWIN  wind  model  as  applied  is  its  inability  to  account  for  frictional  dissi¬ 
pation  of  the  winds  as  a  hurricane  makes  landfall.  Prior  to  reaching  the  land, 
excessive  wind  speeds  create  an  artificially  high  storm  surge  on  the  right-hand 
side  of  the  hurricane  causing  overprediction  by  the  model.  As  winds  rotate 
overland  without  reduction  in  magnitude,  a  diminished  surge  occurs  on  the 
left-hand  side  of  the  hurricane.  Despite  the  fact  that  wind  stress  values  com¬ 
puted  by  the  HURWIN  model  are  generally  too  high  in  the  region  of  hurricane 
landfall,  the  uncalibrated  storm  surge  computations  generated  by  the  ADCIRC 
model  are  quite  reasonable  when  computed  over  the  east  coast  domain. 


Conclusions 

Previously,  little  effort  has  been  made  to  understand  the  relationship 
between  domain  size,  boundary  conditions,  and  the  resulting  physics  for  hurri¬ 
cane  storm  surge  modeling  in  continental  margin  waters.  In  this  domain  size 
sensitivity  study,  storm  surge  computations  over  three  domain  sizes  subject  to 
two  different  open  ocean  boundary  forcings  are  compared  to  detennine  the 
influence  of  domain  size  on  the  computed  response. 

Results  of  this  investigation  clearly  illustrate  that  stomi  surge  model 
domains  which  are  largely  situated  on  the  continental  shelf  and  whose  domain 
size  is  limited  relative  to  the  size  of  the  storm  significantly  underestimate  the 
primary  stoim  surge  response.  Significant  stonu  surge  occurs  in  the  vicinity  of 
open  ocean  boundaries  of  such  domains  when  water  is  pushed  up  on  the  shelf 
by  hurricane  winds.  Therefore,  appropriate  boundary  conditions  are  difficult, 
if  not  impossible,  to  specify  for  these  domains  which  are  small  compared  to 
the  storm  scale  and  which  have  large  expanses  of  cross-shelf  open  ocean  boun¬ 
daries.  Consequently,  a  significantly  underestimated  storm  surge  response 
results  when  using  a  small  domain  such  as  the  Florida  coast  domain.  Despite 
this,  most  storm  surge  modeling  efforts  to  date  have  used  similarly  sized  conti¬ 
nental  shelf  domains,  e.g.  SLOSH  model  applications  by  the  National  Weather 
Service  (Jarvinen  and  Lawrence  1985;  Shaffer,  Jelesnianski,  and  Chen  1986; 
Jelesnianski,  Chen,  and  Shaffer  1992;  Chen,  Shaffer,  and  Kim  1993). 

The  computed  storm  surge  response  over  a  domain  which  encompasses  the 
Gulf  of  Mexico  exhibits  oscillatory  behavior  due  to  the  existence  of  resonant 
modes  or  surge  forerunner  in  the  Gulf  of  Mexico  basin.  These  resonant  modes 
in  the  Gulf  of  Mexico  are  well-documented  and  can  be  easily  excited/ 
influenced  by  a  number  of  factors:  numerical  discretization,  start-up  condi¬ 
tions,  interior  domain  forcing  functions,  and/or  the  boundary  conditions.  This 
study  reinforces  the  findings  of  others  that  resonant  modes  within  the  Gulf  of 
Mexico  are  quite  sensitive  to  boundary  forcing  functions  specified  over  the 
Gulf  of  Mexico  domain.  Consequently,  the  stonn  surge  response  computed 
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over  the  Gulf  of  Mexico  domain  may  not  capture  the  physics  associated  with 
the  hurricane  forerunner. 

The  east  coast  domain,  which  includes  the  western  North  Atlantic  Ocean, 
the  Gulf  of  Mexico,  and  the  Caribbean  Sea,  leads  to  convergent  predictions  of 
both  the  primary  storm  surge  and  surge  forerunner  effect.  A  hurricane  pro¬ 
gresses  through  the  domain  generating  and  propagating  storm  surge  in  a  natu¬ 
ral  and  realistic  fashion.  The  inclusion  of  contiguous  basins  allows  proper 
setup  of  basin  resonant  modes  and  facilitates  the  realistic  propagation  of  storm 
surge  throughout  the  domain  onto  the  continental  shelf,  where  development  of 
storm  surge  is  most  critical.  The  main  advantage  of  the  east  coast  domain  is 
that  the  open  boundaries  lie  within  the  deep  Atlantic  Ocean  and  are  far  from 
the  intricate  processes  occurring  on  the  continental  shelf  and  within  the  Gulf  of 
Mexico  basin  in  response  to  the  storm.  As  was  demonstrated,  sensitivity  of 
the  coastal  response  to  open-ocean  meteorological  boundary  forcing  is  mini¬ 
mal.  Furthermore,  tidal  elevations  appropriate  at  the  deep  Atlantic  Ocean 
boundary  in  the  east  coast  domain  are  readily  extracted  from  global  tidal  mod¬ 
els.  Comparisons  of  storm  surge  elevations  computed  over  the  east  coast 
domain  to  observed  water  levels  verify  that  the  model  response  is  reliable. 
These  comparisons  also  point  out  the  importance  of  accurate  meteorological 
forcing  and  possible  shortcomings  in  the  HURWIN  wind  model. 

Response  characteristics  of  a  storm  surge  model  of  the  continental  margin 
are  profoundly  influenced  by  the  domain  size  and  the  associated  boundary 
conditions.  The  amplitude  of  the  primary  surge  varied  significantly  over  the 
three  domains  examined.  Furthermore,  entire  flow  phenomena  may  be  mere 
model  artifacts  that  appear  and  disappear.  In  fact,  a  significant  Helmholtz 
mode  could  be  easily  excited  in  the  Gulf  of  Mexico  model  domain  while  this 
is  not  the  case  in  either  of  the  other  domains.  Thus,  in  order  to  make 
meaningful  statements  about  physics,  the  modeler  should  be  reasonably  certain 
about  what  level  of  convergence  has  been  achieved  relative  to  domain  size. 
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Discretization  of  the  hydrodynamic  domain  is  just  as  important  as  domain 
size  selection  in  the  formulation  of  a  numerical  storm  surge  model.  The  dis¬ 
crete  form  of  the  hydrodynamic  model  equations  presented  in  Chapter  2  are 
applicable  at  individual  points.  Therefore,  a  discrete  representation  of  the 
computational  domain  is  constructed  by  dividing  the  domain  into  a  mesh  of 
nodd  points  to  form  a  grid.  Computed  stom  surge  elevations  and  velocities 
are  obtained  at  each  nodal  point  within  this  grid  upon  solution  of  the  model 
equations.  The  numerical  solution  procedure  implemented  at  each  grid  node 
utilizes  information  from  surrounding  nodes.  As  a  result,  the  degree  of  spatial 
resolution  provided  by  a  grid  of  discrete  points  significantly  influences  the 
computed  storm  surge  response  of  the  model. 

Too  often,  grids  are  constructed  using  subjective  criteria  and  their  discreti¬ 
zations  are  based  on  computational  constraints.  Only  recently  have  efforts 
been  made  to  determine  if  such  grids  are  convergent  with  respect  to  numerical 
model  computations  of  tidal,  wind-driven,  and  large-scale  baroclinic  circulation 
(Johns  et  al.  1983a,b;  Le  Provost  and  Vincent  1986;  Bennett  and  Campbell 
1987;  Dietrich,  Roache,  and  Marietta  1990;  Lardner  and  Song  1992;  Piacsek 
and  Allard  1993;  Dietrich  1993;  Westerink,  Luettich,  and  Muccino  1994; 
Westerink  et  al.,  in  preparation  (a);  Luettich  and  Westerink  1994). 

In  most  storm  surge  modeling  applications,  though,  no  rigorous  studies  of 
grid  convergence  are  evident,  thou^  grid  resolution  is  recognized  as  being 
important  for  accurate  storm  surge  prediction.  Recall  that  large  domains 
extending  into  the  deep  ocean  are  desirable  for  accurate  storm  surge  prediction. 
Such  domains  require  variably  spaced,  and  often  unstructured,  discretizations 
to  meet  computational  constraints.  Near  coastal  continuous  grid  compression 
(Dube,  Sinha,  and  Roy  1986;  Johns  et  al.  1983a)  is  one  technique  implemented 
to  provide  varying  degrees  of  resolution  throughout  the  model  domain.  How¬ 
ever,  restrictions  on  the  skewness  and  maximum  cell-to-ceU  size  ratios  severely 
limit  the  benefits  of  such  an  approach.  The  use  of  nested  grids  is  a  common 
way  of  handling  resolution  requirements  (e.g.  Flather  (1984);  Dendrou,  Moore, 
and  Myers  (1985);  Al-Rabeh,  Eunay,  and  Cekirge  (1990);  Hubbert,  Leslie,  and 
Manton  (1990)).  In  the  nested  grid  approach,  computations  of  storm  surge 
over  a  large  coarse  grid  are  interpolated  onto  the  boundaries  of  a  local  highly 
refined  grid.  Simulation  of  the  storm  surge  then  continues  separately  over  the 
smaller  grid.  Nested  grids  are  usually  not  coupled  during  a  given  simulation 
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and  nodal  spacing  between  nested  grids  generally  differs  by  a  factor  of  two  to 
three  with  the  minimum  spacing  in  the  finest  grid  ranging  from  0.5  km 
(.3  s  miles)  to  10  km  (6.2  s  miles).  A  coupled,  nested  grid  system  would  more 
closely  approach  the  grid  flexibility  of  the  finite  element  method  but  may  not 
be  as  convenient.  Regardless  of  the  gridding  technique  used,  little  has  been 
done  to  determine  if  the  spatial  resolution  provided  by  the  grid  discretization 
leads  to  accurate  storm  surge  predictions. 

The  focus  of  this  chapter  is  to  systematically  determine  an  optimal  structure 
for  the  domain  discretization  which  leads  to  accurate  hurricane  storm  surge 
computations.  The  grid  structure  sought  is  one  that  provides  only  the  spatial 
resolution  necessary  to  capture  pertinent  storm  surge  physics  and  does  not 
overdiscretize.  This  is  particularly  important  if  model  efficiency  is  to  be  pre¬ 
served  when  using  large  computational  domains  as  suggested  by  the  results  of 
the  domain  size  sensitivity  study.  Findings  from  the  grid  studies  performed 
herein  are  applicable  to  any  of  the  grid  discretization  strategies  currently  in 
use,  though  a  variably  spaced,  unstructured  grid  is  optimal  when  using  very 
large  domains. 

For  this  study  of  grid  structure,  the  placement  of  nodal  points  throughout 
the  domain  is  examined  with  respect  to  the  degree  of  spatial  resolution 
required  for  accurate  representation  of  storm  surge  generation.  The  affect  of 
grid  spacing  and  coastline  resolution  on  hurricane  storm  surge  calculations  is 
also  investigated.  Errors  in  the  primary  surge  computed  over  grids  which  have 
different  nodal  spacings  and  coastline  resolutions  are  related  to  the  grid  spac¬ 
ing  as  well  as  bathymetric  changes  and  characteristics  of  the  hurricane  forcing. 

Investigation  of  grid  discretization  requirements  for  hurricane  storm  surge 
calculations  involves  the  computation  of  storm  surge  elevations  over  an  ideal¬ 
ized  domain  using  various  regular  and  variably  graded  grid  discretizations. 
Meteorological  forcing  is  provided  by  four  synthetic  hurricanes.  Results  of 
this  grid  convergence  study  lead  to  guidelines  for  grid  structure  and  discretiza¬ 
tion  strategies  when  modeling  hurricane  storm  surge  generation. 


Grid  Descriptions 


Domain  and  coastline  definitions 

A  rectangular  domain  with  dimensions  of  2,500  by  3,000  km  (8,200  by 
9,840  miles)  is  constructed  to  correspond  to  the  areal  extent  of  the  east  coast 
domain  used  throughout  Chapter  4.  The  rectangular  domain  defined  does  not 
include  resonant  basins  such  as  the  Gulf  of  Mexico  or  Caribbean  Sea  found 
within  the  east  coast  domain.  Instead  this  domain  more  closely  represents  the 
Atlantic  Ocean  off  the  eastern  shore  of  the  United  States.  A  land  boundary 
representing  the  coastline  lies  along  the  left-most  3,000-km  (9,840-mile)  length 
of  the  rectangular  domain  and  open  ocean  conditions  apply  at  the  remaining 
domain  boundaries.  A  representative  bathymetry  profile,  shown  in  Figure  40 
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Continental  Shelf 


DISTANCE  FROM  COAST  (km) 


Figure  40.  Bathymetry  profile  for  the  rectangular  domain  used  in  the  study  of  grid  structure 
relative  to  storm  surge  predictions 


relative  to  the  rectangular  domain,  is  synthesized  from  ocean  depth  cross  sec¬ 
tions  recorded  off  the  coasts  of  Virginia  and  Florida  in  the  western  North 
Atlantic  Ocean  as  well  as  average  slopes  and  distances  suggested  in  the  litera¬ 
ture  for  portions  of  the  bathymetry  profile  (Weyl  1970;  Weisberg  and  Parish 
1974;  Pickard  and  Emery  1982;  Bearmad  1989;  Thurman  1974). 

Two  different  representations  of  the  coastline  are  considered.  The  first  is  a 
simple  straight-line  coast.  This  particular  land  boundary  eliminates  the  effect 
of  coastline  variability  on  storm  surge  generation.  A  second  coastline  is  con¬ 
structed  using  a  sine  wave  variation,  which  has  an  amplitude  of  30  km 
(18.6  statute  miles)  and  a  wavelength  of  100  km  (62.1  s  miles)  along  the 
length  of  the  shore  and  is  given  by: 
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where  (Xg,  yg)  are  the  coordinates,  in  kilometers,  of  the  extreme  northwest 
land  boundary  point  in  the  domain,  x  is  the  horizontal  coordinate  of  a  point  at 
the  coast,  and  y  is  the  distance,  in  kilometers,  southward  along  the  shoreline. 
The  sinusoidal  deviations  of  the  coastline  are  intended  to  portray  in  an  ideal¬ 
ized  way  the  highly  irregular  profile  of  the  actual  shoreline.  For  simulations 
over  domains  having  this  sinusoidal  land  boundary,  the  importance  of  spatial 
coastline  variability  on  storm  surge  generation  is  examined. 

Over  domains  having  a  sinusoidal  coastline,  the  bathymetry  profile  taken 
from  the  land  boundary  is  represented  by  the  horizontal  cross  section  shown  in 
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Figure  40.  This  results  in  a  north-south  sinusoidal  bathymetric  variation  due  to 
the  deviation  of  coastal  node  positions.  In  the  ocean,  the  influence  of  shore¬ 
line  bathymetric  variations  does  not  extend  beyond  the  continental  shelf.  Thus, 
the  north-south  sinusoidal  variation  of  bathymetry  is  attenuated  over  a  distance 
equal  to  three  times  the  amplitude  of  the  sine  wave  variation  at  the  coast 
(~  90  km  (55.9  miles)).  This  ensures  that  effects  of  coastline  variation  on  the 
bathymetry  profile  remain  on  the  continental  shelf. 


Discretizations 

Six  grid  discretizations  of  the  rectangular  domain  form  the  basis  of  this  grid 
convergence  study.  Three  uniform  grids,  designated  as  GOl,  G02,  and  G03 
and  shown  in  Figure  41  have  regular  nodal  spacings  of  50  km  (31.1  statute 
miles),  25  km  (15.5  statute  mUes,  and  12.5  km  (7.8  statute  miles),  respectively. 
Three  variably  graded  grids,  VGOl,  VG02,  and  VG03  (shown  in  Figure  42) 
have  nodal  spacings  ranging  between  12.5  km  (7.8  miles)  and  50  km 
(31.1  miles).  Grid  VGOl  (Figure  42a)  has  the  minimum  nodal  spacing  of 
12.5  km  (7.8  miles)  confined  to  a  region  near  the  coastline.  The  mesh  spacing 
then  rapidly  increases  to  50  km  (31.1  miles)  at  a  point  approximately  50  km 
(31.1  miles)  offshore.  In  grid  VG02  (Figure  42b)  refinement  at  the  12.5-km 
(7.8-mile)  level  extends  from  the  shoreline  to  25  km  (15.5  miles)  offshore  and 
25-km  (15.5-mile)  spacing  covers  the  remainder  of  the  continental  shelf  with 
50-km  (31.1-mile)  spacing  elsewhere.  The  final  graded  grid,  VG03  (shown  in 
Figure  42c),  has  extensive  uniform  resolution  of  12.5  km  (7.8  miles)  over  the 
entire  continental  shelf,  shelf  break,  and  portions  of  the  continental  slope  out  to 
a  depth  of  approximately  1,000  m  (3,280  ft).  Spacing  over  the  remainder  of 
the  deep  ocean  increases  rapidly  from  25  km  to  50  km  (15.5  to  31.1  miles). 

Understanding  the  influence  of  coastline  resolution  is  an  important  part  of 
the  study  of  domain  discretization.  Three  coastline  representations  named  Cl, 
C2,  and  C3  and  shown  in  Figure  43  have  5-,  9-,  and  17-point-per-wavelength 
resolutions  of  the  sinusoidally  varying  coastline,  respectively.  These  three 
coastline  discretizations  are  merged  with  the  six  grid  discretizations  previously 
described  to  form  a  total  of  fourteen  grids.  Table  4  summarizes  the  distin¬ 
guishing  features  of  the  discretization  associated  with  each  of  the  fourteen 
grids. 


Error  Analysis 


Grid  comparisons 

Four  series  of  grid  comparisons  are  defined  using  the  fourteen  grids  out¬ 
lined  in  Table  4.  Each  of  these  comparisons  is  designed  to  examine  a  particu¬ 
lar  aspect  of  the  grid  discretization  in  relation  to  its  impact  on  storm  surge 
predictions.  Computations  of  storm  surge  over  each  grid  within  a  series  are 
compared  to  the  storm  surge  elevations  computed  over  a  “truth”  grid 
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Figure  41 .  Enlargements  of  the  (a)  50-km,  (b)  25-km,  and  (c)  12.5-km  uniform  discretizations 
for  grids  G01 ,  G02,  and  G03,  respectively 
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Figure  42.  Enlargements  of  the  variably  graded  discretizations  for  grids  (a)  VG01 ,  (b)  VG02, 
and  (c)  VG03,  respectively 
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Figure  43.  Sinusoidally  varying  coastlines  (a)  Cl,  (b)  C2,  and  (c)  C3,  having 
resolutions  of  5,  9,  and  17  points  per  wavelength,  respectively 


Table  4 

Characteristics  of  Rectangular  Grid  Discretizations 

Name 

Domain 

Coastline 

Number  of  Nodes 

Structure 

Spacing 

Type 

Refinement 

G01 

Regular 

50  km 

Straight 

.... 

3,111 

G02 

Regular 

25  km 

Straight 

.... 

12,221 

G03 

Regular 

1 2.5  km 

Straight 

.... 

48,441 

VGOI 

Graded 

Variable 

Straight 

.... 

3,717 

VG02 

Graded 

Variable 

Straight 

.... 

4,014 

VG03 

Graded 

Variable 

Straight 

— 

5,639 

G01_C1 

Regular 

50  km 

Sinusoidal 

5  pts 

3,171 

G02_C1 

Regular 

25  km 

Sinusoidal 

5  pts 

12,401 

G03_C1 

Regular 

12,5  km 

Sinusoidal 

5  pts 

48,920 

G02_C2 

Regular 

25  km 

Sinusoidal 

9  pts 

12,401 

G03_C3 

Regular 

12.5  km 

Sinusoidal 

17  pts 

49,040 

VG01_C3 

Graded 

Variable 

Sinusoidal 

17  pts 

4,136 

VG02_C3 

Graded 

Variable 

Sinusoidal 

17  pts 

4,374 

VG03_C3 

Graded 

Variable 

Sinusoidal 

17  pts 

5,910 
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designated  within  that  series.  The  first  series  of  grid  comparisons  investigates 
grid  spacing  effects  without  the  influence  of  a  spatially  variable  coastline. 
Thus,  grids  GOl,  G02,  G03,  VGOl,  VG02,  and  VG03,  all  of  which  have 
straight  coastlines,  make  up  the  Series  1  grid  comparisons.  The  “truth”  grid 
solution  for  Series  1  is  obtained  using  grid  G03.  A  second  series  of  grids. 
Series  2,  which  compares  grids  G01_C1  and  G02_C1  to  a  “truth”  grid  solution 
over  grid  G03_C1,  demonstrates  the  effect  of  a  spatially  varying  shoreline. 
These  grids  all  have  uniform  mesh  spacings  and  a  5-point-per-wavelength  rep¬ 
resentation  of  the  sinusoidally  varying  coastline.  Another  set  of  grids. 

Series  3,  compares  grid  G03_C1  to  grid  G03_C3  and  grid  G02_C1  to 
grid  G02_C2,  in  order  to  assess  the  consequence  of  coastline  refinement  alone. 
A  final  series  of  comparisons,  named  Series  4,  between  grids  G01_C1, 
G02_C2,  VG01_C3,  VG02_C3,  and  VG03_C3  with  computations  over 
grid  G03_C3  designated  as  a  “truth”  solution,  combines  the  influence  of  grid 
spacing  and  coastline  resolution  on  storm  surge  computation.  Table  5  summa¬ 
rizes  the  grid  comparisons  just  outlined  and  notes  the  aspect  of  the  grid  discre¬ 
tization  which  is  isolated  by  each  grid  comparison  grouping.  The  truth  grid 
designated  within  each  series  of  grid  comparisons  is  listed  in  boldface  in 
Table  5. 


Table  5 

Definitions  of  Grid  Comparisons 

Series 

Grids  Compared 

Objective 

1 

G01  G02  VG01  VG02  VG03  G03 

Grid  spacing 

2 

G01_C1  G02_C1  G03_C1 

Coastline  variation 

3 

G03_C1  G03_C3;  G02_C1  G02_C2 

Coastline  refinement 

4 

G01_C1  G02_C2  VG01_C3 

VG02_C3  VG03_C3  G03_C3 

Combination 

Storm  surge  computations 

For  the  storm  surge  simulations  pertinent  to  this  grid  convergence  study, 
four  synthetic  hurricanes  (HOll,  H012,  H013,  and  H031)  described  in  detail  in 
Chapter  3  serve  as  the  meteorological  forcing  over  all  fourteen  grids.  Impor¬ 
tant  properties  of  these  synthetic  hurricanes  are  repeated  in  Table  6  for  con¬ 
venience.  In  the  simulations  of  storm  surge,  wind  stress  and  pressure  forcing 
are  applied  on  the  interior  of  the  domain  and  tidal  forcing  is  neglected  both  on 
the  interior  and  at  the  open  ocean  boundaries.  Additionally,  an  inverted 
barometer  condition  is  specified  at  the  open  ocean  boundaries  of  the  domain. 

For  this  series  of  simulations,  only  the  finite  amplitude  terms  are  not 
included  in  the  governing  equations  due  to  instabilities  caused  by  near  drying 
elements.  For  consistency,  model  parameters  are  the  same  as  those  used  in  the 
domain  comparison  study  (i.e.,  the  bottom  friction  coefficient  is  constant  and 
equal  to  0.003  and  the  GWCE  parameter  Tq  is  defined  equal  to  0.001).  In 
addition,  model  parameters  are  identical  for  aU  simulations  in  the  grid 
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Table  6 

Characteristics  of  Synthetic  Hurricanes 

Hurricane 

Name 

(km) 

Velocity 

(km/hr) 

Hurricane  Path 
(Relative  to 
the  coast) 

Length  of 
Hurricane 
(days) 

Time  of 

Landfall 

(hrs) 

H011 

30 

15 

Perpendicular 

8.75 

192 

H012 

60 

15 

Perpendicular 

8.75 

192 

H013 

30 

25 

Perpendicular 

6.00 

126 

H031 

30 

25 

Parallel 

20.50 

— 

convergence  study  so  that  comparisons  between  the  computed  storm  surge  for 
different  grid  discretizations  are  possible. 

Simulations  are  spun  up  from  homogeneous  initial  conditions  using  a  1-day 
ramp  in  time.  An  identical  ramp  function  of  1  day  length  is  applied  to  the 
wind  and  pressure  forcing  as  well  as  the  inverted  barometer  boundary  condi¬ 
tion.  Actual  simulation  periods  run  6  hr  longer  than  the  length  of  each  hurri¬ 
cane  listed  in  Table  6.  This  is  due  to  a  6-hr  ramp-up  period  at  the  start  of 
simulation  during  which  initial  hurricane  wind  and  pressure  forcings  are  held 
stationary.  Thereafter,  storm  surge  computations  use  the  time-varying  wind 
and  pressure  fields.  A  time-step  of  45  sec  is  used  throughout  the  simulation 
period.  No  calibration  or  tuning  of  parameters  is  performed  in  either  the 
weather  model  or  in  the  hydrodynamic  model. 


Error  measures 

The  objective  of  the  error  analysis  is  to  determine  the  magnitude  and  loca¬ 
tion  of  errors  in  the  predicted  storm  surge  throughout  the  computational 
domain.  Using  this  information,  the  level  of  resolution  required  to  minimize 
these  errors  is  determined.  In  general,  this  error  analysis  guides  the  placement 
of  nodes  throughout  the  domain  such  that  the  discretization  of  the  domain  uses 
a  minimum  number  of  points,  yet  yields  an  accurate  representation  of  the 
storm  surge  generation  in  the  coastal  ocean. 

Storm  surge  predictions  computed  over  the  most  finely  uniform  discretized 
grid  within  each  series  of  grid  comparisons,  shown  in  bold  in  Table  5,  are  con¬ 
sidered  “truth”  solutions  and  are  used  to  assess  errors  in  the  storm  surge  com¬ 
putations  performed  over  the  remaining  grids  within  a  series. 

Maximum  overprediction  (positive)  and  underprediction  (negative)  errors  in 
the  storm  surge  are  computed  over  each  grid  as  a  function  of  time.  The  abso¬ 
lute  error  is  calculated  simply  as  the  difference  between  computed  storm  surge 
values  and  the  “true”  storm  surge  elevations  at  each  point  in  time.  The 
relative  error  normalizes  the  absolute  error  with  respect  to  the  highest  surge 
elevation  over  the  entire  domain  at  the  specific  point  in  time  being  considered. 
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An  inclusion  of  ocean  bathymetry  at  each  point  in  time  considered  indicates 
the  location  of  the  maximum  over-  and  underprediction  and  helps  identify 
where  grid  resolution  is  lacking. 


Discussion  of  Resuits 

For  the  discussion  of  storm  surge  prediction  errors  that  follows,  knowledge 
of  the  storm  surge  elevation  profile  for  each  hurricane  is  helpful  in  understand¬ 
ing  the  error  patterns  that  develop.  Figures  44a  -  44d  depict  the  maximum 
storm  surge  generated  over  the  domain  at  a  specific  time  and  the  depth  at  the 
location  of  maximum  water  elevation  for  meteorological  forcing  by  Hurri¬ 
canes  HOll,  H012,  H013,  and  H031.  The  surge  profiles  in  Figures  44a  -  44d 
are  taken  from  computations  over  grid  VG03  and  will  be  referenced  in  discuss¬ 
ing  the  storm  surge  prediction  errors  presented  for  each  series  of  grid  compari¬ 
sons.  Note  the  oscillatory  behavior  of  the  surge  elevations  in  deep  water. 

This  is  an  artifact  of  the  insufficient  resolution  provided  over  the  deep  ocean 
and  will  be  discussed  further  in  subsequent  sections. 


Series  1  Comparisons:  G01,  G02,  VG01,  VG02,  VG03 

For  Hurricane  HOll,  maximum  overprediction  and  underprediction  errors 
computed  over  uniform  grids  GOl  and  G02  are  shown  in  Figures  45a  and  45b, 
respectively.  Up  to  about  174  hr,  while  the  storm  is  predominately  located  in 
the  deep  ocean  (Figure  44a),  the  oscillating  error  pattern  seen  in  Figures  45a 
and  45b  is  associated  with  an  insufficient  resolution  in  deep  water  to  capture 
the  inverted  barometer  which  forms  (i.e.,  the  pressure  forcing  function  for 
Hurricane  HOll,  which  has  a  spatial  scale  of  30  km  (18.6  s  miles).  This 
oscilUating  error  pattern  is  also  insufficiently  resolved  over  grids  having  uni¬ 
form  spacings  of  50  km  (31.1  statute  miles)  and  25  km  (15.5  statute  miles), 
and  results  in  an  aliased  response  which  fluctuates  depending  on  the  position 
of  the  hurricane  eye  relative  to  the  nearest  node).  This  variably  aliased  behav¬ 
ior  is  also  apparent  in  the  storm  surge  elevation  profile  shown  in  Figure  44a 
while  Hurricane  HOI  1  is  in  the  deep  ocean. 

As  the  storm  moves  onto  the  continental  shelf  at  about  180  hr,  the  absolute 
errors  increase  sigruficantly.  This  error  increase  corresponds  to  the  dramatic 
increase  in  storm  surge  as  winds  push  water  up  onto  the  shelf  against  the 
coast.  Also  note  the  large  increase  in  the  relative  error  which  coincides  with  a 
decreasing  absolute  error  and  a  declining  peak  surge  as  the  storm  makes  land¬ 
fall.  This  behavior  of  the  relative  error  indicates  that  processes  which  dissipate 
the  peak  surge  are  much  more  rapid  than  those  which  dissipate  the  errors  in 
the  storm  surge.  A  comparison  of  the  error  magnitudes  between  Figures  45a 
and  45b  demonstrates  that  halving  the  nodal  spacing  leads  to  a  reduction  of  the 
error  in  deep  waters  by  a  factor  of  two  and  a  threefold  reduction  of  the  errors 
generated  on  the  continental  shelf.  Peak  errors  over  grids  GOl  and  G02  are 
concentrated  on  the  continental  shelf  nearest  the  shoreline  suggesting  that  a 


Chapter  5  Influence  of  Grid  Structure 


Figure  44.  Maximum  storm  surge  elevations  over  the  domain  and  the  depth  at  the  location  of 
maximum  surge  computed  over  grid  VG03  for  (a)  Hurricane  H01 1  and  (b)  Hurri¬ 
cane  H012  forcing  (Continued) 
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Figure  44.  (Concluded) 
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Figure  45.  Maximum  over-  and  underprediction  errors  of  the  storm  surge  generated  by  Hurri¬ 
cane  H01 1  as  computed  over  grids  (a)  G01 ,  (b)  G02,  (c)  VG01 ,  (d)  VG02,  and 
(e)  VG03  (Sheet  1  of  3) 
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Figure  45.  (Sheet  2  of  3) 
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variably  graded  grid  may  provide  an  optimal  grid  structure  for  minimization  of 
the  storm  surge  prediction  errors. 

For  the  three  variably  graded  grids  (VGOl,  VG02,  and  VG03)  maximum 
storm  surge  prediction  errors  and  depths  at  the  location  of  the  maximum  errors 
for  these  grids  are  presented  in  Figures  45c  -  45e  relative  to  Hurricane  HO  11 
forcing.  For  each  of  the  graded  grids,  errors  in  the  storm  surge  computations 
over  the  deep  ocean  are  the  same  as  those  computed  over  grid  GOl  due  to  an 
identical  grid  spacing  of  50  km.  As  the  storm  moves  onto  the  continental 
shelf  at  about  180  hr,  the  largest  oveiprediction  errors  computed  relative  to 
grid  VGOl,  seen  in  Figure  45c,  initially  occur  over  the  continental  slope,  a 
shift  from  the  predominate  shoreline  errors  computed  over  grids  GOl  and  G02 
at  the  same  times.  The  additional  resolution  provided  at  the  coastline  in 
grid  VGOl  improves  storm  surge  predictions  at  the  coast  until  the  time  of  peak 
surge  and  hurricane  landfall  at  192  hr.  Following  the  passage  of  Hurri¬ 
cane  HOll  onto  land,  maximum  errors  over  grid  VGOl  remain  on  the  conti¬ 
nental  shelf  but  are  removed  from  the  vicinity  of  the  coastline.  In  comparison 
to  grid  G02,  maximum  errors  over  grid  VGOl  are  reduced  by  a  factor  of  four 
due  to  the  additional  refinement  provided  at  the  coastline.  In  general,  the 
maximum  errors  computed  over  grid  VGOl  are  relatively  uniform  in  time  and 
space  and  are  significantly  less  than  those  computed  over  either  grid  GOl 
or  G02. 
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Maximum  error  levels  for  graded  grid  VG02,  presented  in  Figure  45d,  indi¬ 
cate  that  the  extreme  errors  are  located  on  the  continental  shelf  away  from  the 
coastline  except  at  the  time  of  peak  surge  when  maximum  errors  remain  at  the 
shore.  Increased  shelf  resolution  over  grid  VG02  reduces  error  levels  below 
that  of  grid  VGOl  but  not  significantly.  Over  VG03,  which  has  extensive 
resolution  over  the  continental  shelf,  shelf  break,  and  slope,  error  reduction 
extends  into  deeper  waters.  Results  from  these  error  analyses  indicate  that 
storm  surge  calculations  in  nonresonant  basins,  where  surge  foreranners  are  not 
developed,  do  not  require  the  same  high  level  of  grid  resolution  over  the  shelf 
break  and  slope  as  tidal  computations  (Westerink,  Luettich,  and  Muccino 
1994;  Westerink  et  al.,  in  preparation  (b);  Westerink,  Luettich,  and  Hagen 
1994). 

Maximum  over-  and  underprediction  errors  in  the  computed  storm  surge 
generated  by  forcing  from  Hurricane  H012  are  shown  in  Figures  46a  -  46e  for 
each  grid  in  Series  1.  Recall  that  Hurricane  H012  has  a  spatial  scale  double 
that  of  Hurricane  HOll  (i.e.,  is  60  km  (32.4  n.m.)).  Notice  while  the 
storm  is  in  the  deep  ocean,  at  times  less  than  174  hr,  the  variably  aliased  error 
pattern  is  still  evident  over  the  uniform  grid  GOl  but  the  magnitude  of  the 
errors  is  half  of  those  associated  with  Hurricane  HOll  forcing  during  the  same 
time  period.  Over  grid  G02,  maximum  errors  computed  in  the  deep  ocean  are 
negligible  as  seen  in  Figure  46b.  Reduced  deepwater  errors  computed  over 
grids  GOl  and  G02  are  attributed  to  the  larger  spatial  scale  of  Hurricane  HOI 2 
(i.e.,  more  nodal  points  are  available  to  represent  the  features  of  a  larger  scale 
storm  resulting  in  an  increased  resolution  of  the  inverted  barometer  effect  as 
well  as  wind  forcing).  For  grid  G02,  the  uniform  spacing  of  25  km  is  approxi¬ 
mately  one  half  the  spatial  scale  of  Hurricane  H012. 

Maximum  overprediction  errors  over  grid  GOl  occur  at  the  coastline  long 
before  the  storm  reaches  the  continental  shelf.  As  the  storm  progresses  onto 
the  continental  shelf  at  about  180  hr,  the  magnitude  of  the  absolute  errors  over 
grid  GOl  is  somewhat  higher  than  those  computed  for  Hurricane  HOll  forcing, 
a  consequence  of  the  increased  peak  surge  generated  by  the  large-scale  storm. 
However,  when  mesh  spacing  is  halved  as  for  grid  G02,  the  same  error  trends 
are  exhibited  for  Hurricane  H012  forcing  as  were  seen  for  forcing  by  Hurri¬ 
cane  HOll.  Namely,  a  threefold  reduction  in  error  over  the  continental  shelf  is 
observed  in  comparing  computed  errors  over  grids  GOl  and  G02  shown  in 
Figures  46a  and  46b.  Maximum  relative  errors  associated  with  Hurri¬ 
cane  H012  forcing  computed  over  all  five  grids  in  Series  1  are  the  same  as  or 
slightly  less  than  those  recorded  for  Hurricane  HOI  1  forcing.  Generally,  a 
large-scale  storm  has  more  nodal  points  available  than  a  smaller  scale  storm  to 
resolve  the  same  features  of  the  hurricane  forcing.  So  naturally,  for  a  given 
mesh  spacing,  more  accurate  predictions  result  for  the  larger  storm  event. 

Relative  to  Hurricane  H012  forcing,  errors  computed  over  grids  VGOl, 
VG02,  and  VG03  while  the  storm  is  in  the  deep  ocean  are  the  same  as  those 
computed  over  grid  GOl  due  to  the  identical  grid  spacings  of  50  km  as  seen  in 
Figures  46c  -  46e.  Over  grid  VGOl,  maximum  error  values,  shown  in  Fig¬ 
ure  46c,  remain  in  deeper  waters  until  12  hr  prior  to  the  time  of  landfall  when 
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Figure  46.  Maximum  over-  and  underprediction  errors  of  the  storm  surge  generated  by  Hurri¬ 
cane  H012  as  computed  over  grids  (a)  G01,  (b)  GG2,  (c)  VG01,  (d)  VG02,  and 
(e)  VG03  (Sheet  1  of  3) 
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maximum  prediction  errors  shift  to  the  continental  slope.  At  the  time  of  hurri¬ 
cane  landfall,  192  hr,  maximum  errors  move  to  the  coastline.  Unlike  errors 
associated  with  grids  GOl  and  G02,  the  increased  refinement  near  the  coastline 
in  grid  VGOl  eliminates  sizable  errors  at  the  coastline  while  the  storm  is  over 
deeper  waters.  The  resolution  provided  over  graded  grid  VGOl  reduces  errors 
to  an  acceptable  uniform  level  throughout  the  domain  as  seen  in  Figure  46c. 
Further  discretization  over  the  shelf  leads  to  little  change  in  the  magnitude  or 
location  of  prediction  errors  as  seen  in  Figure  46d  for  grid  VG02.  Similar  to 
the  results  presented  for  Hurricane  HOll  forcing,  resolution  over  the  entire 
shelf  and  beyond  the  shelf  break  shifts  maximum  errors  back  into  the  deep 
ocean  as  seen  in  Figure  46e  for  grid  VG03.  Error  levels  over  grid  VG03  dur¬ 
ing  the  24-hr  period  surrounding  the  time  of  peak  surge  are  reduced  to  mag¬ 
nitudes  seen  over  deeper  water. 

Application  of  the  meteorological  forcing  from  Hurricane  HO  13  over  the 
Series  1  grid  discretizations  demonstrates  the  effect  of  an  increased  forward 
speed  on  grid  resolution  requirements.  The  forward  speed  of  Hurricane  HO  13 
is  approximately  1.5  times  that  of  Hurricane  HOll.  The  storm  surge  hydro¬ 
graph  for  Hurricane  HO  13  in  Figure  44c  exhibits  a  peak  surge  which  is  much 
larger  and  occurs  over  a  significantly  shorter  time  period  than  the  surge 
produced  by  Hurricane  HOll.  As  a  result  of  the  increased  storm  surge  genera¬ 
tion,  the  maximum  over-  and  underprediction  errors  computed  over  each  grid 
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Figure  47.  Maximum  over-and  underprediction  errors  of  the  storm  surge  generated  by  Hurri¬ 
cane  H013  as  computed  over  grids  (a)  G01,  (b)  G02,  (c)  VG01,  (d)  VG02,  and 
(e)  VG03  (Sheet  1  of  3) 
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in  Series  1,  shown  in  Figures  47a  -  47e,  deviate  slightly  upwards  from  the 
errors  associated  with  Hurricane  HO  11  over  the  same  grids. 

Since  the  spatial  scale  of  Hurricane  HOI  3  is  the  same  as  that  for  Hurri¬ 
cane  HOll,  errors  occurring  while  the  storm  is  in  the  deep  ocean  are  similar  to 
those  reported  for  Hurricane  HOll  and  indicate  under-resolution  of  the  inverted 
barometer  effect  in  deep  water.  Maximum  storm  surge  prediction  errors  over 
grids  GOl  and  G02  occur  at  the  shoreline  18  hr  in  advance  of  the  storm  entry 
onto  the  continental  shelf  at  128  hr.  A  reduction  in  grid  spacing  by  a  factor  of 
2  halves  the  errors  in  the  deep  ocean  while  on  the  continental  shelf  errors  are 
diminished  by  a  factor  of  2  to  3  when  comparing  maximum  errors  relative  to 
grids  GOl  and  G02  shown  in  Figures  47a  and  47b.  Further  reduction  of  the 
nodal  spacing  at  the  coastline  in  grid  VGOl  leads  to  errors  which  are  dimin¬ 
ished  by  more  than  a  factor  of  3  from  those  over  grid  G02  and  the  locations  of 
these  errors  coincide  with  the  path  of  the  storm  as  it  moves  onto  the  continen¬ 
tal  shelf.  Only  minor  differences  are  evident  in  a  comparison  of  the  magni¬ 
tudes  and  locations  of  the  maximum  errors  shown  in  Figures  47c  and  47d  for 
grids  VGOl  and  VG02.  The  extension  of  maximum  errors  into  deeper  waters 
and  the  substantial  reduction  in  the  magnitude  of  the  maximum  error  over 
grid  VG03,  shown  in  Figure  47e,  indicate  that  errors  over  the  continental  slope 
may  be  more  significant  for  the  case  of  forcing  by  a  relatively  fast-moving 
hurricane. 
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Figure  48.  Maximum  over-  and  underprediction  errors  of  the  storm  surge  generated  by  Hurri¬ 
cane  H031  as  computed  over  grids  (a)  G01,  (b)  G02,  (c)  VG01,  (d)  VG02,  and 
(e)  VG03  (Sheet  1  of  3) 
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Maximum  errors  computed  over  the  grids  in  Series  1  for  hurricane  forcing 
which  remains  parallel  to  the  coastline,  Hurricane  H031,  are  shown  in  Fig¬ 
ures  48a  -  48e.  The  maximum  surge  profile  for  Hurricane  H031,  depicted  in 
Figure  44d,  shows  more  sustained  but  reduced  surge  elevations.  While  the 
hurricane  is  in  the  deep  ocean,  errors  computed  over  all  grids  are  similar  to 
those  recorded  for  forcing  by  Hurricane  HOll.  As  seen  in  Figures  48a 
and  48b,  which  compare  storm  surge  prediction  errors  over  grids  GOl 
and  G02,  the  maximum  error  in  deep  waters  decreases  by  a  factor  of  2  when 
grid  spacing  is  halved.  As  the  storm  approaches  the  continental  shelf,  maxi¬ 
mum  errors  over  grids  GOl  and  G02  are  concentrated  at  the  shoreline  and 
remain  there  over  24  hr  after  the  storm  has  returned  to  deeper  waters.  Gen¬ 
erally,  on  the  continental  shelf  absolute  errors  computed  over  the  uniform 
grids  GOl  and  G02  relative  to  H031  forcing  are  less  than  those  experienced 
with  respect  to  Hurricane  HOI  1  forcing.  This  observation  is  correlated  to  the 
decreased  storm  surge  heights  generated  by  Hurricane  H031.  Both  the  abso¬ 
lute  and  relative  errors  reduce  threefold  when  nodal  spacing  is  halved  between 
grids  GOl  and  G02. 

Graded  grid  VGOl,  which  has  twice  the  resolution  of  grid  G02  at  the  coast¬ 
line,  yields  maximum  prediction  errors  which  are  located  primarily  offshore  as 
the  storm  moves  parallel  to  the  coastline  as  seen  in  Figure  48c.  These  errors 
have  magnitudes  which  are  one  half  the  value  of  those  computed  over 
grid  G02.  While  the  storm  is  on  the  continental  shelf,  errors  over  grid  VG02, 


90 


Chapter  5  Influence  of  Grid  Structure 


shown  in  Figure  48d,  are  essentially  uniform  and  error  magnitudes  are  reduced 
slightly  from  those  associated  with  grid  VGOl.  Additionally,  the  offshore 
locations  of  the  maximum  error  are  quite  similar  between  grids  VGOl 
and  VG02.  Upon  finely  resolving  the  entire  shelf  and  slope,  storm  surge  pre¬ 
diction  errors  are  further  reduced  and  their  locations  are  shifted  to  the  deep 
ocean  as  seen  in  Figure  48e  for  grid  VG03.  Generally,  though,  the  level  of 
resolution  provided  by  grid  VG02,  which  produces  low  uniform  errors  over  the 
continental  shelf,  is  satisfactory. 

The  maximum  values  of  the  computed  errors  over  both  the  deep  waters  and 
continental  shelf  regions  of  the  ocean  are  presented  in  Tables  7  -  10  for  each 
grid  in  Series  1  relative  to  meteorological  forcing  from  Hurricanes  HOll, 

H012,  H013,  and  H031.  The  maximum  error  values  recorded  in  Tables  7  -  10 
relative  to  each  hurricane  forcing  clearly  exhibit  the  trends  discussed  previ¬ 
ously  when  grid  resolution  is  increased  by  a  factor  of  two  among  grids  GOl, 
G02,  and  VGOl.  Over  deep  waters,  a  reduction  of  the  grid  spacing  by  half 
diminishes  the  error  in  the  computed  storm  surge  by  a  factor  of  two.  On  the 
continental  shelf  and  in  coastal  regions,  the  error  reduction  is  generally  three¬ 
fold.  Maximum  storm  surge  prediction  errors  computed  over  the  variably 
graded  grids  are  minimal  and  essentially  uniform  throughout  the  domain. 
Furthermore,  differences  in  the  errors  associated  with  each  of  the  graded 
grids  VGOl,  VG02,  and  VG03  are  often  slight  with  regard  to  the  error 
magnitudes. 


Table  7 

Maximum  Errors  Over  Grid  Series  1  for  Hurricane  H01 1 

Grid 

Overprediction 

Underprediction 

Absolute  (cm) 

Relative 

Absolute  (cm) 

Relative 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

G01 

4-9 

153 

10-18 

59 

8-22 

68 

12-51 

61 

G02 

1-4 

57 

3-8 

28 

8-12 

30 

12-22 

24 

VG01 

4-9 

11 

10-18 

10 

8-22 

4 

12-51 

8 

VG02 

4-9 

9 

10-18 

5 

8-22 

2 

12-51 

4 

VG03 

4-9 

- 

10-18 

- 

8-22 

- 

12-51 

- 

In  comparing  the  errors  among  Tables  7-10  for  different  hurricane  forc¬ 
ings,  the  relative  improvement  of  the  deep-ocean  resolution  for  forcing  by  the 
more  expansive  storm.  Hurricane  H012,  is  evident.  In  Tables  8  and  9,  the 
increased  values  of  the  maximum  error  associated  with  Hurricanes  HO  12 
and  HOI 3  are  correlated  to  the  significant  storm  surge  generated  by  these  hur¬ 
ricanes.  Relative  errors  computed  over  each  grid  in  Tables  7-10  are  strik¬ 
ingly  similar  regardless  of  the  hurricane  forcing  applied.  This  behavior  of  the 
relative  error  reinforces  that  resolution  at  the  coastline,  not  the  character  of  the 
hurricane  forcing,  is  the  predominate  factor  affecting  the  accurate  prediction  of 
hurricane  storm  surge. 
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Table  8 

Maximum  Errors  Over  Grid  Series  1  for  Hurricane  H012 

Grid 

Overprediction 

Underprediction 

Absolute  (cm) 

Relative 

Absolute  (cm) 

Relative 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

G01 

2-4 

179 

6-8 

52 

8-11 

94 

12-21 

122 

G02 

2 

55 

4 

21 

2-4 

35 

4-8 

15 

VG01 

2-4 

10 

6-8 

5 

8-11 

5 

12-21 

3 

VG02 

2-4 

9 

6-8 

5 

8-11 

5 

12-21 

3 

VG03 

2-4 

- 

6-8 

- 

8-11 

- 

12-21 

- 

Table  9 

Maximum  Errors  Over  Grid  Series  1  for  Hurricane  HOI 3 

Grid 

Overprediction 

Underprediction 

Absolute  (cm) 

Relative 

Absolute  (cm) 

Relative 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

G01 

5-10 

167 

9-16 

72 

7-22 

91 

12-58 

56 

G02 

2 

52 

4 

37 

4-9 

41 

10-12 

28 

VG01 

5-10 

11 

9-16 

9 

7-22 

9 

12-58 

60 

VG02 

5-10 

9 

9-16 

8 

7-22 

9 

12-58 

6 

VG03 

5-10 

“ 

9-16 

- 

7-22 

- 

12-58 

- 

Table  10 

Maximum  Errors  Over  Grid  Series  1  for  Hurricane  H031 

Grid 

Overprediction 

Underprediction 

Absolute  (cm) 

Relative 

Absolute  (cm) 

Relative 

Deep 

Shelf 

1 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

G01 

4-8 

95 

9-15 

65 

23 

57 

58 

40 

G02 

2-4 

32 

4-8 

31 

8-12 

27 

12-25 

28 

VG01 

4-8 

32 

9-15 

20 

23 

12 

69 

12 

VG02 

4-8 

9 

9-15 

15 

23 

29 

69 

9 

VG03 

4-8 

- 

9-15 

- 

23 

- 

69 

- 
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Series  2  Comparisons:  G01_C1,  G02_C1 


The  second  series  of  grid  comparisons  examines  the  influence  of  a  spatially 
varying  shoreline.  The  two  grids  studied,  G01_C1  and  G02_C1,  have  uniform 
spacings  of  50  km  (31.1  statute  miles)  and  25  km  (15.5  statute  miles),  respec¬ 
tively,  and  use  a  5-point-per- wavelength  representation  of  the  sinusoidal  coast¬ 
line.  Maximum  overprediction  and  underprediction  errors  computed  relative  to 
forcing  by  Hurricane  HOI  1  are  shown  in  Figures  49a  and  49b.  Errors  com¬ 
puted  while  the  storm  is  in  the  deep  ocean  are  similar  to  those  computed  over 
grids  using  a  straight-line  coast  indicating  that  the  influence  of  the  coastline 
does  not  extend  past  the  continental  shelf.  At  156  hr,  prior  to  the  storm  entry 
onto  the  shelf,  the  overprediction  errors  computed  with  respect  to  grid  G01_C1 
(Figure  49a)  are  already  located  at  the  coastline.  Absolute  error  values  com¬ 
puted  over  both  grid  G01_C1  and  grid  G02_C1  are  nearly  twice  the  value  of 
those  errors  computed  over  grids  having  the  same  domain  discretization  but 
without  spatial  coastline  variations;  i.e.,  grids  GOl  and  G02.  The  sinusoidally 
varying  coastline  magnifies  the  storm  surge  elevations  generated  at  the  shore. 
In  comparing  grid  G01_C1  to  grid  G02_C1,  a  reduction  in  the  grid  spacing  by 
a  factor  of  2  leads  to  a  fourfold  decrease  in  the  absolute  errors  at  the  coastline. 
Relative  errors  are  similar  in  magnirnde  to  those  computed  over  grids  GOl 
and  G02  in  Series  1.  Table  11  records  the  maximum  overprediction  and 
underprediction  error  values  for  the  two  grids  in  Series  2  subject  to  forcing  by 
Hurricane  HOll. 


Series  3  Comparisons:  G03_C1,  G02_C1 

The  impact  of  coastline  resolution  is  examined  by  considering  the  storm 
surge  errors  produced  when  comparing  a  5-point  to  a  17-point-per-wavelength 
resolution  of  the  sinusoidally  varying  coastline  (i.e.,  a  comparison  between 
grids  G03_C1  and  G03_C3)  and  by  comparing  a  5-point  to  a  9-point-per- 
wavelength  resolution  of  the  sinusoidally  varying  coastline  (i.e.,  a  comparison 
between  grids  G02_C1  and  G02_C2).  The  applied  meteorological  forcing  is 
from  Hurricane  HOll.  Since  each  of  the  grids  compared  has  identical  and  fine 
discretizations  everywhere  except  at  the  coastline,  errors  throughout  the 
domain  are  negligible,  as  seen  in  Figures  50a  and  50b.  At  the  coastline,  how¬ 
ever,  over-  and  underprediction  errors  have  extreme  magnitudes  ranging  from 
20  cm  to  40  cm,  (0.6  ft  to  1.3  ft),  as  seen  in  Table  12,  which  summarizes  the 
maximum  error  values  for  the  Series  3  grid  comparisons.  The  discrepancy 
between  storm  surge  predictions  made  using  different  coastline  representations 
clearly  illustrates  the  importance  of  coastline  resolution  in  the  accurate  predic¬ 
tion  of  hurricane  storm  surge. 


Series  4  Comparisons:  G01_C1,  G02_C2,  VG01_C3,  VG02_C3, 
VG03_C3 

The  Series  4  grid  comparisons  include  the  influences  of  grid  spacing,  spa¬ 
tial  coastline  variation,  and  coastline  resolution  on  storm  surge  prediction.  The 
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TIME  (hrs) 


(a) 


TIME  (hrs) 


(b) 


Figure  49.  Maximum  over-  and  underprediction  errors  of  the  storm  surge  generated  by  Hurri¬ 
cane  H011  as  computed  over  grids  (a)  G01_C1,  and  (b)  G02_C1 
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Table  11 

Maximum  Errors  Over  Grid  Series  2  for  Hurricane  H011 

Grid 

Overprediction 

Underprediction 

Absolute  (cm) 

Relative 

Absolute  (cm) 

Relative 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

G01_C1 

4-9 

320  , 

10-18 

53 

8-22 

188 

12-51 

71 

G02_C1 

1-4 

80 

3-8 

21 

8-12 

41 

12-22 

19 

storm  surge  generated  by  each  of  the  four  synthetic  hurricanes  is  simulated 
over  all  grids  in  Series  4.  The  maximum  errors  in  computed  storm  surge  are 
presented  for  Hurricane  HOll  forcing  in  Figures  51a  -  51e,  for  Hurri¬ 
cane  H012  forcing  in  Figures  52a  -  52e,  for  Hurricane  HOI 3  forcing  in  Fig¬ 
ures  53a  -  53e,  and  for  Hurricane  H031  forcing  in  Figures  54a  -  54e.  Analysis 
of  the  overprediction  and  underprediction  errors  reconfirms  the  general  obser¬ 
vations  noted  in  the  grid  comparisons  for  Series  1-3.  Namely,  while  the  storm 
is  in  the  deep  ocean,  prediction  errors  exhibit  a  variably  aliased  behavior  indi¬ 
cating  under-  resolution  of  the  inverted  barometer  forcing  as  well  as  the  wind 
forcing  function.  For  the  larger  scale  storm.  Hurricane  HO  12,  the  relative 
deep-ocean  resolution  is  essentially  halved  in  comparison  to  the  size  of  the 
storm  and  consequently,  errors  are  reduced  twofold  as  seen  in  Figures  52a 
and  52b. 

A  sinusoidally  varying  coastline  amplifies  the  absolute  errors  computed 
over  the  uniform  grids  G01_C1  and  G02_C2  to  at  least  twice  the  value  of  the 
errors  computed  over  the  straight-line  coast  grids  GOl  and  G02  (see  Fig¬ 
ures  51a  and  51b:  52a  and  52b;  53a  and  53b;  54a  and  54b).  As  grid  spacing 
over  the  uniform  grids  is  cut  by  a  factor  of  2,  the  absolute  and  relative  errors 
in  computed  storm  surge  are  reduced  twofold  in  the  deep  ocean  and  as  much 
as  three  times  over  the  continental  shelf.  Relative  maximum  errors  remain 
similar  in  magnitude  to  those  reported  in  the  Series  1  grid  comparisons  as  seen 
in  Tables  13-16  which  present  the  maximum  error  values  for  each  grid  in 
Series  4  subject  to  forcing  from  the  four  synthetic  hurricanes.  Generally, 
though,  the  magnitude  of  the  errors  for  the  Series  4  grid  comparisons  are 
somewhat  elevated  over  the  error  magnitudes  recorded  for  the  Series  1 
comparisons. 

The  effect  of  spatial  variation  at  the  shoreline  is  an  increased  generation 
and  accumulation  of  water  along  the  coast.  The  variably  graded  discretization 
of  grid  VG01_C3,  which  includes  a  highly  resolved  coastline,  leads  to  a 
decrease  in  the  magnitude  of  the  prediction  errors  over  the  continental  shelf  by 
a  factor  of  2  to  3  as  seen  in  Tables  13-16.  This  reduction  in  error  is  due  to 
the  combination  of  a  finely  discretized  coastline  and  the  doubling  of  resolution 
near  the  shore  to  12.5  km  (see  Figures  51c,  52c,  53c,  and  54c).  Despite  the 
diminished  prediction  errors,  error  magnitudes  over  grid  VG01_C3  remain  at 
unacceptable  levels  between  20  and  30  cm.  Over  graded 
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Figure  50.  Maximum  over-  and  underprediction  errors  of  the  storm  surge  generated  by  Hurri¬ 
cane  H011  as  computed  over  grids  (a)  G03_C1,  and  (b)  G02_C1  when  compared 
to  grids  G03_C3  and  G02_C2,  respectively 
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Table  12 

Maximum  Errors  Over  Grid  Series  3  for  Hurricane  H011 


Overprediction 

Underprediction 

Grid 

Absoiute  (cm) 

Relative 

Absolute  (cm) 

Relative 

Deep 

Sheif 

Deep 

Sheif 

Deep 

Shelf 

Deep 

Shelf 

G03_C1 

~ 

38 

- 

11 

- 

21 

- 

12 

G02_C1 

- 

18 

- 

8 

-- 

30 

- 

10 

grids  VG02_C3  and  VG03_C3,  maximum  prediction  errors  are  further  reduced 
from  the  levels  computed  over  VG01_C3,  but  error  magnitudes  at  the  coastline 
remain  significant.  While  the  storms  are  on  the  continental  shelf,  computed 
errors  over  grid  VG03_C3  persist  at  the  coastline  unlike  the  deepwater  loca¬ 
tions  of  the  prediction  errors  computed  relative  to  the  straight  coastline  graded 
grid,  VG03.  Upon  resolution  of  the  shelf,  shelf  break,  and  portions  of  the  con¬ 
tinental  slope,  it  has  been  previously  shown  that  errors  over  grid  VG03  shifted 
to  deep  waters,  whereas  errors  over  grid  VG03_C3  remain  predominately  at 
the  coastline.  For  each  of  the  hurricane  forcings  applied  in  the  Series  4  grid 
comparisons,  storm  surge  prediction  errors  persist  at  the  coastline  regardless  of 
the  dscretization  strategy.  Thus,  the  grid  comparisons  in  Series  4  conclusively 
demonstrate  that  the  spatial  variation  of  the  coastline  dominates  the  generation 
of  storm  surge.  Resolution  of  and  around  the  shoreline  is  critical  for  accurate 
prediction  of  hurricane  storm  surge. 


Conclusions 

A  study  of  the  influence  of  grid  structure  on  hurricane  storm  surge  predic¬ 
tion  demonstrates  that  discretization  of  the  computational  domain  can  signifi¬ 
cantly  affect  the  storm  surge  elevation  computed  in  the  coastal  region.  A 
decrease  in  mesh  spacing  by  a  factor  of  2  halves  the  errors  over  the  deep 
ocean  and  results  in  a  two-  to  threefold  reduction  of  errors  in  the  coastal 
region  for  the  storms  studied.  This  relationship  between  errors  in  the  com¬ 
puted  storm  surge  and  mesh  spacing  is  observed  until  error  magnitudes  stabi¬ 
lize  and  are  relatively  uniform  throughout  the  domain.  At  this  point,  storm 
surge  predictions  over  a  particular  grid  discretization  are  considered  acceptable 
within  a  defined  error  limit.  Over  the  deep  ocean,  errors  in  the  predicted 
storm  surge  are  minimized  by  using  a  nodal  spacing  which  is  approximately 
one  half  the  spatial  scale  of  the  hurricane.  In  coastal  areas,  significant  refine¬ 
ment  both  of  the  coastline  geometry  and  nearshore  regions  is  necessary  to 
eliminate  excessive  error  in  storm  surge  predictions.  Determination  of  the 
exact  grid  spacing  required  in  coastal  areas  depends  primarily  on  the  complex¬ 
ity  of  the  coastline  detail  and  the  extent  of  very  shallow  waters. 
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Figure  51.  Maximum  over-  and  underprediction  errors  of  the  storm  surge  generated  by  Hurri¬ 
cane  H011  as  computed  over  grids  (a)  G01_C1,  (b)  G02_C2,  (c)  VG01_C3, 

(d)  VG02_C3,  and  (e)  VG03_C3  (Sheet  1  of  3) 
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(C) 


(d) 


Figure  51 .  (Sheet  2  of  3) 
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Figure  51. 


(Sheet  3  of  3) 


High  levels  of  resolution  over  the  continental  shelf  break  and  slope  do  not 
appear  to  be  necessary  for  the  storm  surge  calculations  performed  in  this 
nonresonant  basin  grid  structure  study.  Rather,  the  finest  level  of  resolution  is 
needed  adjacent  to  the  coast  in  shallow  waters.  Furthermore,  refinement  of  the 
coastline  detail  is  shown  to  be  very  important  for  accurate  prediction  of  storm 
surge  elevations.  Significant  errors  are  recorded  when  coastline  detail  alone  is 
under-resolved. 

In  relation  to  the  parameters  of  hurricane  forcing  (i.e.,  path,  forward 
velocity,  and  spatial  scale)  resolution  near  the  coastline  is  still  the  most  signifi¬ 
cant  factor  in  storm  surge  computations.  Faster  moving  and/or  larger  storms 
magnify  the  errors  in  under-resolved  regions  of  the  domain,  particularly  near 
the  coast.  One  possible  exception  is  for  stonns  which  track  parallel  to  the 
coastline.  In  addition  to  the  generation  of  stonn  surge  at  the  coast,  these 
storms  produce  significant  surge  elevations  over  the  entire  continental  shelf 
and  slope.  As  a  result,  higher  resolution  of  the  continental  shelf  break  and 
slope  may  improve  storm  surge  predictions  for  hurricanes  approaching  along  a 
path  parallel  to  the  shoreline.  The  level  of  resolution  required  over  the  conti¬ 
nental  shelf  and  slope  will  depend  on  the  areal  extent  of  the  storm  and  the 
relative  position  of  the  path  to  the  shelf  break. 
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Figure  52.  Maximum  over-  and  underprediction  errors  of  the  storm  surge  generated  by  Hurri¬ 
cane  H012  as  computed  over  grids  (a)  G01_C1.  (b)  G02_C2,  (c)  VG01_C3, 

(d)  VG02_C3,  and  (e)  VG03_C3  (Sheet  1  of  3) 
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Figure  52.  (Sheet  2  of  3) 
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Figure  52.  (Sheet  3  of  3) 


The  disparity  between  grid  resolution  requirements  in  the  deep  ocean  and 
coastal  regions  suggests  that  a  variably  graded  grid  structure  is  most  appropri¬ 
ate.  For  the  grid  study  and  error  analyses  conducted,  a  variable  grid  structure 
which  has  extensive  refinement  in  the  nearshore  region  leads  to  low  uniform 
errors  throughout  the  domain.  Since  discretization  errors  can  never  be  com¬ 
pletely  eliminated,  uniform  errors  over  a  grid  are  desired.  Even  more  impor¬ 
tantly,  errors  in  the  prediction  of  storm  surge  generation  near  the  land  and  on 
the  continental  shelf  must  be  minimized.  An  appropriate  grid  discretization  is 
constructed  by  considering  the  phenomenon  of  interest  (i.e.,  primary  surge 
and/or  resonant  mode  behavior),  the  accuracy  constraints  for  storm  surge  pre¬ 
dictions  in  the  region  of  interest,  desired  error  levels  throughout  the  remainder 
of  the  domain,  and  computational  effort. 

In  addition  to  providing  low  uniform  errors  over  the  domain,  the  variable 
grid  structure  requires  minimum  computational  effort  which  is  particularly 
noticeable  when  implemented  in  conjunction  with  large  domains.  The  flexibil¬ 
ity  of  the  finite  element  approach  used  to  discretize  the  grids  within  this  study 
leads  to  easy  incorporation  of  coastline  detail  and  nodal  densities  which  range 
over  an  order  of  magnitude.  The  variation  in  nodal  density  arises  from  the 
significant  refinement  provided  in  shallow  coastal  areas  where  storm  surge 
generation  is  important  and  in  regions  of  complex  coastline  detail  and/or  bath¬ 
ymetric  change,  in  conjunction  with  the  coarse  discretizations  over  the  deep 
ocean  where  processes  occur  more  gradually  and  are  of  less  interest. 
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Figure  53.  Maximum  over-  and  underprediction  errors  of  the  storm  surge  generated  by  Hurri¬ 
cane  H013  as  computed  over  grids  (a)  G01_C1,  (b)  G02_C2,  (c)  VG01_C3, 

(d)  VG02_C3,  and  (e)  VG03_C3  (Sheet  1  of  3) 
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Figure  53.  (Sheet  2  of  3) 
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Figure  53.  (Sheet  3  of  3) 


Implementation  of  a  graded  grid  discretization  over  large  domains,  coupled 
with  the  efficiency  of  the  finite  element  method,  leads  to  a  discrete  problem 
that  remains  well  within  computational  limits. 

Furthermore,  the  error  analysis  procedure  implemented  to  assess  the  quality 
of  predicted  storm  surge  elevations  offers  a  straightforward  means  to  evaluate 
the  performance  of  a  particular  grid  discretization.  The  magnitude  and  loca¬ 
tion  of  errors  in  the  computed  stonn  surge  are  readily  obtained  and  can  assist 
in  the  construction  of  an  optimal  grid  structure  which  yields  accurate  storm 
surge  predictions. 
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(a) 


(b) 


Figure  54.  Maximum  over-  and  underprediction  errors  of  the  storm  surge  generated  by  Hurri¬ 
cane  H031  as  computed  over  grids  (a)  G01_C1,  (b)  G02_C2,  (c)  VG01_C3, 

(d)  VG02_C3,  and  (e)  VG03_C3  (Sheet  1  of  3) 
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Figure  54.  (Sheet  3  of  3) 


Table  13 

Maximum  Errors  Over  Grid  Series  4  for  Hurricane  H011 

Grid 

Overprediction 

Underpredlction 

Absolute  (cm) 

Relative 

Absolute  (cm) 

Relative 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

G01_C1 

4 

-9 

331 

10-18 

53 

8-22 

191 

12-51 

80 

G02_C2 

1-4 

80 

3-8 

20 

8  12 

48 

12-22 

20 

VG01_C3 

4 

-9 

32 

10-18 

21 

8-22 

20 

12-51 

11 

VG02_C3 

4 

-9 

19 

10-18 

7 

8-22 

17 

12-51 

11 

VG03_C3 

4 

-9 

19 

10-18 

7 

8-22 

14 

12-51 

6 
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Table  14 

Maximum  Errors  Over  Grid  Series  4  for  Hurricane  HOI 2 

Grid 

Overprediction 

Underprediction 

Absoiute  (cm) 

Relative 

Absolute  (cm) 

Relative 

Deep 

Sheif 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

G01_C1 

2-4 

413 

6-8 

48 

8-11 

310 

12-21 

60 

G02_C2 

2 

92 

4 

21 

2-4 

63 

4-8 

21 

VG01_C3 

2-4 

34 

6-8 

12 

8-11 

19 

12-21 

7 

VG02_C3 

2-4 

17 

6-8 

6 

8-11 

21 

12-21 

10 

VG03_C3 

2-4 

20 

6-8 

8 

8-11 

16 

12-21 

7 

Table  15 

Maximum  Errors  Over  Grid  Series  4  for  Hurricane  HOI 3 

Grid 

Overprediction 

Underprediction 

Absolute  (cm) 

Relative 

Absolute  (cm) 

Relative 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

G01_C1 

5-10 

346 

9-16 

68 

7-22 

227 

12-58 

71 

G02_C2 

2 

81 

4 

22 

4-9 

89 

10-12 

40 

VG01_C3 

5-10 

40 

9-16 

16 

7-22 

21 

12-58 

8 

VG02_C3 

5-10 

28 

9-16 

13 

7-22 

21 

12-58 

18 

VG03_C3 

5-10 

21 

9-16 

11 

7-22 

20 

12-58 

10 

Table  16 

Maximum  Errors  Over  Grid  Series  4  for  Hurricane  H031 

Grid 

Overprediction 

Underprediction 

Absolute  (cm) 

Relative 

Absolute  (cm) 

Relative 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

Deep 

Shelf 

G01_C1 

4-8 

205 

9-15 

75 

23 

246 

58 

109 

G02_C2 

2-4 

62 

4-8 

29 

8-12 

90 

12-25 

53 

VG01_C3 

4-8 

34 

9-15 

21 

23 

31 

69 

20 

VG02_C3 

4-8 

25 

9-15 

15 

23 

23 

69 

13 

VG03_C3 

4-8 

22 

9-15 

12 

23 

21 

69 

17 
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6  Application  of  Domain  Size 
and  Gridding  Strategy 


The  domain  size  sensitivity  study  and  the  investigation  of  grid  resolution 
effects  have  established  several  guidelines  for  selection  of  a  domain  and  con¬ 
struction  of  the  domain  discretization  associated  with  a  numerical  storm  surge 
model. 

To  properly  represent  the  generation  and  propagation  of  storm  surge  into 
coastal  areas,  a  very  large  domain  is  required.  This  domain  includes  the 
coastal  region  of  interest,  contiguous  basins,  and  extends  out  of  resonant  basins 
into  the  deep  open  ocean.  A  domain  of  this  size  allows  proper  setup  of  reso¬ 
nant  modes/surge  forerunner.  Furthermore,  specification  of  elevations  at  the 
open  ocean  boundary  is  simplified  and  has  little  perceptible  influence  on  the 
model  response  in  coastal  regions. 

An  optimal  discretization  of  the  domain  assumes  a  graded  structure  with  a 
high  level  of  resolution  provided  in  near  coastal  regions.  Nodal  spacing  in  the 
deep  ocean  is  less  dense  and  should  approximate  one  half  the  spatial  scale  of 
the  hurricane  forcing.  Significant  coastline  detail  must  be  included  in  the  dis¬ 
cretization  as  well  if  accurate  storm  surge  predictions  are  to  result. 

An  application  of  the  numerical  modeling  strategy  just  outlined  with  respect 
to  domain  size  and  grid  structure  is  presented.  Forcing  from  two  historical 
hurricanes  is  applied  over  a  large  domain  having  a  variably  graded  grid  discre¬ 
tization.  Errors  in  the  storm  surge  predicted  over  this  grid  discretization  are 
identified  by  comparison  to  a  “truth”  grid  solution.  The  “truth”  grid  for  these 
applications  no  longer  has  uniform  mesh  spacing  but  is  a  variably  graded 
discretization,  defined  by  a  fourfold  refinement  of  the  resolution  provided  in 
the  original  study  grid.  Differences  between  storm  surge  elevations  computed 
over  the  study  grid  and  the  “truth”  grid  are  examined  to  determine  the  magni¬ 
tude  and  location  of  the  maximum  storm  surge  prediction  errors.  Spatial  dis¬ 
tributions  of  prediction  errors  are  also  utilized  to  assess  the  performance  of  the 
study  grid  in  accurately  representing  the  storm  surge  response  in  coastal 
regions. 

The  applications  presented  demonstrate  the  utility  of  the  large  domain  and 
graded  grid  structure  advocated  for  use  in  conjunction  with  numerical  storm 
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surge  models.  Furthermore,  the  procedure  applied  to  assess  the  level  of  grid 
convergence  with  respect  to  storm  surge  computations  is  easily  implemented 
and  can  be  used  effectively  to  assist  in  the  construction  of  an  optimal  grid 
discretization. 


Domain  and  Grid  Construction 

Two  historical  hurricanes,  Kate  (November  1985)  and  CamiUe  (August 
1969),  will  serve  as  the  meteorological  forcing  for  these  storm  surge  model 
applications.  Details  regarding  the  paths  and  characteristics  of  Hurricanes  Kate 
and  Camille  are  specified  in  Chapter  3.  Each  of  these  hurricanes  made  landfall 
at  points  along  the  Gulf  of  Mexico  coast  in  the  United  States  (i.e..  Hurricane 
Kate  at  Panama  City,  FL,  and  Hurricane  CamiUe  near  Biloxi,  MS). 

An  appropriate  domain,  then,  is  one  that  centers  on  regions  surrounding 
Biloxi,  MS,  and  Panama  City,  FL,  but  also  includes  the  Gulf  of  Mexico  basin 
and  extends  into  the  deep  regions  of  the  western  North  Atlantic  Ocean.  The 
domain  selected  is  based  on  the  east  coast  domain,  which  encompasses  the 
western  North  Atlantic  Ocean,  the  Caribbean  Sea,  and  the  Gulf  of  Mexico  as 
described  and  used  in  Chapter  4.  The  inclusion  of  greater  coastline  detail  in 
the  region  of  Biloxi,  MS,  shown  in  Figure  55,  (e.g.  the  addition  of  Lake 
Pontchartrain,  Louisiana  west  of  Biloxi,  MS)  is  the  only  change  made  to  the 
east  coast  domain.  The  deep  Atlantic  Ocean  boundary  remains  along  the 
60  °W  meridian  and  aU  other  boundaries  are  defined  by  the  eastern  coastlines 
of  North,  Central,  and  South  America.  Topography  within  the  domain  is  the 
same  as  that  of  the  east  coast  domain  except  in  the  region  surrounding  Biloxi, 
MS.  In  this  area,  bathymetry  has  been  updated  using  the  National  Oceanic  and 
Atmospheric  Administration  Digital  U.S.  Coastal  Hydrography  sounding  data¬ 
base.  RecaU  from  the  description  of  the  east  coast  domain  that  the  continental 
shelf  has  depths  ranging  from  an  imposed  minimum  between  3  m  (10  ft)  and 
7  m  (23  ft)  to  130  m  (426  ft)  at  the  shelf  break,  the  continental  slope  has  a 
typical  depth  range  of  130  m  to  3,000  m  (426  ft  to  9,840  ft),  and  depths  of  the 
continental  rise  and  deep  ocean  increase  upwards  from  3,000  m  (9,840  ft)  to 
nearly  8,000  m  (26,240  ft). 

Having  defined  the  model  domain,  a  discrete  representation  of  the  domain, 
grid  SGOl,  is  constructed  using  23,566  nodes  and  43,238  elements.  The  vari¬ 
ably  graded  structure  of  the  SGOl  grid  discretization  yields  nodal  spacings  that 
range  from  a  maximum  of  98  km  (60.9  statute  miles)  in  the  deep  ocean  to  a 
minimum  of  0.5  km  (0.3  statute  miles)  near  coastal  areas.  A  map  depicting 
the  distribution  of  grid  resolution  provided  throughout  the  domain  of  grid 
SGOl  is  shown  in  Figure  56.  As  demonstrated  in  the  study  of  grid  structure  in 
Chapter  5,  capturing  the  inverted  barometer  effect  of  the  hurricane  over  deep 
waters  requires  a  level  of  resolution  which  is  at  least  one  half  the  spatial  scale 
of  the  hurricane.  Throughout  the  period  of  simulation.  Hurricane  Camille 
maintains  a  large  spatial  scale  of  approximately  80  km  (49.7  statute  miles). 
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Figure  55.  Region  of  coastline  detail  surrounding  Biloxi,  MS,  which  has  been  added  to  the 
east  coast  domain  used  throughout  Chapter  4 


The  tracking  of  Hurricane  Camille  begins  near  Cuba  and  follows  a  path  over 
the  Gulf  of  Mexico  where  the  maximum  grid  spacing  is  approximately  40  km 
(24.8  statute  miles).  Thus,  over  deep  waters,  the  resolution  provided  for  Hurri¬ 
cane  Camille  is  indeed  half  the  spatial  scale  of  the  hurricane. 

For  Hurricane  Kate,  the  spatial  scale  of  the  storm  over  the  deep  ocean  is 
highly  variable  and  ranges  from  20  km  to  80  km  (12.4  statute  miles  to 
49.7  statute  miles).  During  the  initial  18  hr  of  Hurricane  Kate  over  the  deep 
western  North  Atlantic  Ocean,  grid  spacing  in  the  vicinity  of  the  path  of  Hurri¬ 
cane  Kate  is  between  50  km  (31.1  statute  miles)  and  60  km 
(37.3  statute  miles).  Consequently,  during  this  period,  the  inverted  barometer 
effect  associated  with  the  pressure  forcing  of  Hurricane  Kate  is  under-resolved. 
However,  as  Hurricane  Kate  moves  towards  Cuba,  grid  spacing  rapidly 
decreases  to  30  km  (18.6  statute  miles)  or  less  and  at  entrances  to  the  Gulf  of 
Mexico  and  Caribbean  Sea  (which  are  important  for  resonant  mode  set-up),  the 
grid  spacing  drops  to  15  km  (9.3  statute  miles)  and  below.  Except  for  very 
deep  portions  of  the  western  North  Atlantic  Ocean,  the  spatial  scale  of  Hurri¬ 
cane  Kate  is  generally  well  resolved  using  the  discretization  provided  by  grid 
SGOl. 

Another  consideration  in  the  construction  of  a  grid  is  the  representation  of 
the  coastline.  Coastline  detail  plays  an  important  role  in  the  generation  of 
storm  surge.  Results  from  the  grid  sensitivity  study  in  Chapter  5  clearly  dem¬ 
onstrate  that  representation  of  the  shoreline  contour  is  critical  to  achieving 
accurate  computations  of  hurricane  storm  surge.  Highly  detailed  shoreline 
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Figure  56.  Distribution  of  nodal  spacing,  in  kilometers,  over  grid  SG01 


coordinates  obtained  from  the  CIA  database  and  NOAA  bathymetric  charts  are 
used  to  define  the  shoreline  profile  of  grid  SGOl  along  the  gulf  coast.  The 
representation  of  coastline  detail  and  the  discretization  associated  with 
grid  SGOl  in  the  gulf  coast  region  near  the  vicinity  of  hurricane  landfall 
locations  are  shown  in  Figure  57.  A  map  of  grid  spacing  for  this  same  region 
is  presented  in  Figure  58.  In  the  near-shore  regions  of  grid  SGOl,  grid  spacing 
ranges  from  0.5  km  (0.3  statute  miles)  to  2  km  (1.2  statute  miles).  The  finest 
resolution  in  Figure  58  is  located  at  the  entrances  to  inlets,  where  setup  of 
exchange  processes  with  the  coastal  ocean  is  important.  Over  a  majority  of  the 
continental  shelf,  the  grid  resolution  is  less  than  5  km  (3.1  statute  miles). 
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Figure  57,  Coastline  detail  and  discretization  for  gird  SG01  in  the  northeast  Gulf  of  Mexico 


For  the  grid  convergence  analysis  which  follows,  a  highly  refined  “truth” 
grid  (CGOl)  is  also  constructed.  The  discretization  of  grid  CGOl  is  comprised 
of  90,435  nodes  and  172,952  elements,  exactly  a  fourfold  increase  of  the  grid 
resolution  provided  over  grid  SGOl.  Storm  surge  predictions  over  grid  SGOl 
will  be  compared  to  computations  over  grid  CGOl  to  determine  the  relative 
accuracy  of  storm  surge  elevations  computed  over  grid  SGOl. 


Storm  Surge  Simulations 

A  series  of  simulations  were  conducted  to  investigate  the  performance  of 
the  selected  domain  size  and  the  grid  discretization  strategy  implemented  for 
the  prediction  of  hurricane  storm  surge.  Computations  were  made  over 
grids  SGOl  and  CGOl,  previously  discussed.  Hurricane  CamiUe  and  Hurricane 
Kate  are  used  separately  as  the  wind  and  pressure  forcing  over  the  entire 
model  domain.  Along  the  open  ocean  boundary,  an  inverted  barometer  pres¬ 
sure  forcing  is  applied.  Simulations  performed  using  Hurricane  Kate  as  forc¬ 
ing  implement  a  fully  nonlinear  model  formulation  with  the  exception  of  the 
finite  amplitude  terms.  These  terms  are  not  included  in  the  simulations 
because  of  instabilities  caused  by  near-drying  elements  in  coastal  areas.  For 
the  case  of  Hurricane  Camille,  the  model  formulation  eliminates  all  nonlineari¬ 
ties  except  the  bottom  friction  terms.  The  magnitude  of  the  surge  generated  by 
Hurricane  Camille  in  conjunction  with  instabilities  in  the  convective  terms 


Chapter  6  Application  of  Domain  Size  and  Gridding  Strategy 


115 


Figure  58.  Distribution  of  nodal  spacing  in  kilometers  over  grid  SG01  in  the  northeast  Gulf  of 
Mexico 


computed  over  grid  CGOl  have  led  to  this  formulation.  As  in  prior  simula¬ 
tions,  model  parameters  are  specified  such  that  the  bottom  friction  coefficient 
is  constant  and  equal  to  0.003  and  the  GWCE  parameter  Tq  is  defined  equal 
to  0.(X)1. 

Simulations  are  spun  up  from  homogeneous  initial  conditions  using  a  1-day 
ramp  in  time.  Application  of  the  hyperbolic  ramp  function  reduces  the  excita¬ 
tion  of  nonphysical  short  wavelength  frequencies.  An  identical  ramp  function 
of  1  day  length  is  applied  to  the  wind  and  pressure  forcing  as  well  as  the 
inverted  barometer  boundary  condition.  Actual  simulations  for  Hurricane  Kate 
begin  at  12:00  GMT  15  November  1985  and  mn  over  8.25  days  (including  a 
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1-day  ramp-up  period).  For  Hurricane  Camille,  model  computations  begin  at 
12:00  GMT  14  August  1969  and  run  for  8.00  days  (including  the  1-day 
ramp-up  period).  During  the  first  6  hr  of  each  simulation,  the  initial  hurricane 
wind  and  pressure  forcings  are  held  stationary.  Thereafter,  storm  surge 
computations  use  the  time -varying  wind  and  pressure  fields.  The  time-step  for 
simulations  over  grid  SGOl  is  45  sec.  This  time-step  length  is  halved  to 
22.5  sec  for  simulations  using  grid  CGOl  so  that  the  Courant  number  based  on 
wave  celerity  is  identical  to  that  over  grid  SGOl  and  remains  below  1.5  for 
accuracy  considerations  (Westerink  et  al.  1992b).  All  parameters  in  both  the 
hydrodynamic  and  wind  models  remain  uncalibrated. 

Seventy-five  elevation  stations,  shown  in  Figures  59  -  61,  are  placed  along 
the  eastern  and  gulf  coasts  of  the  United  States,  on  the  continental  shelf,  and  at 
the  shelf  break  in  the  northeastern  portion  of  the  Gulf  of  Mexico  near  the 
landfall  regions  of  Hurricanes  Kate  and  Camille.  Comparisons  between  storm 
surge  elevation  solutions  computed  over  grids  SGOl  and  CGOl  will  be  made  at 
these  stations. 

Discossion  of  Convergence  Analysis 


Hurricane  Camille 

An  evaluation  of  the  storm  surge  predictions  computed  over  grid  SGOl 
begins  with  an  examination  of  the  maximum  storm  surge  errors  over  the 
domain  at  specific  points  in  time.  The  location  and  magnitude  of  the  max¬ 
imum  overprediction  and  underprediction  errors  are  shown  in  Figure  62. 
Absolute  and  relative  error  measures  are  the  same  as  those  prescribed  in  Chap¬ 
ter  5  with  the  “tmth”  solution  now  defined  by  storm  surge  elevations  computed 
over  grid  CGOl. 

As  seen  in  Figure  62,  all  maximum  prediction  errors  occur  at  the  coastline 
throughout  the  period  of  simulation.  Furthermore,  the  maximum  errors  exhibit 
peaks  in  magnitude  at  several  times.  The  first  peak  of  the  absolute  overpredic¬ 
tion  error,  approximately  110  cm,  is  seen  at  84  hr,  the  time  of  peak  surge  and 
hurricane  landfall.  The  second  noticeable  increase  in  the  absolute  overpre¬ 
diction  error  occurs  after  138  hr  as  Hurricane  Camille  moves  off  land  and 
returns  to  the  coastal  waters  near  Chesapeake  Bay-Bridge,  Virginia.  Areas 
along  the  shoreline  of  the  eastern  United  States  are  not  under  consideration  and 
thus  have  not  been  discretized  to  the  same  levels  found  in  the  Gulf  of  Mexico. 
The  increase  in  absolute  error  after  138  hr  is  a  direct  consequence  of  this 
under-resolution  along  the  eastern  coast  of  the  United  States. 

Figure  62  exhibits  relative  errors  that  remain  fairly  uniform  during  the 
approach  and  landfall  of  Hurricane  Camille  within  the  Gulf  of  Mexico  and  on 
the  Mississippi  shelf  The  extreme  relative  errors  coincide  with  hurricane 
positions  in  deep  waters  or  over  land.  At  a  time  of  126  hr,  the  relative  error 
maximum  is  due  to  a  combination  of  residual  errors  at  the  coastline  following 
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Figure  60.  Additional  storm  surge  elevation  stations,  shown  as  solid  squares, 
in  and  around  the  Gulf  of  Mexico 


Figure  61 .  Seven  storm  surge  elevation  stations,  shown  as  solid  squares, 
along  the  eastern  United  States  coastline 
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Figure  62,  Maximum  over-  and  underprediction  errors  in  the  storm  surge 
generated  by  Hurricane  Camille  as  computed  over  grid  SG01 


landfall  of  the  hurricane  and  the  presence  of  minimal  surge  elevations  as  the 
hurricane  moves  inland.  Relative  errors  are  at  a  minimum  as  Hurricane 
Camille  approaches  the  coast  and  makes  landfall.  The  unifonn,  low-magnitude 
profile  of  the  relative  prediction  error  for  the  simulation  period  corresponding 
to  Hurricane  Camille's  movement  through  the  Gulf  of  Mexico  indicates  the 
appropriateness  of  the  variably  graded  discretization  in  grid  SGOl.  However, 
high  magnitudes  of  the  absolute  error  at  the  time  of  pealc  surge  warrant  exami¬ 
nation  of  the  spatial  distribution  of  storm  surge  prediction  errors  throughout 
the  domain. 

A  series  of  snapshots  taken  at  different  times  during  the  simulation  and 
shown  in  Figure  63  depict  the  absolute  overprediction  (shown  as  shades  of  red) 
and  underprediction  (shown  in  shades  of  blue)  errors  of  the  storm  surge  calcu¬ 
lated  over  grid  SGOl  as  compared  to  computations  made  over  grid  CGOl.  As 
Hurricane  Camille  approaches  land,  errors  increase  over  an  18-hr  period  sur¬ 
rounding  the  time  of  landfall  but  are  confined  to  a  limited  region  at  the  coast 
near  Biloxi,  MS.  After  84  hr  on  18  August  1969,  errors  build  around  an 
island  directly  in  the  path  of  Hurricane  Camille  as  seen  in  Figure  63c.  An 
enlargement  of  the  discretization  in  this  region,  presented  in  Figure  64,  shows 
relatively  coarse  resolution  around  the  island,  which  would  indeed  cause  larger 
Storm  surge  prediction  errors.  Six  hours  later,  overprediction  of  the  storm 
surge  reaches  the  1-m  level  in  a  localized  area  between  the  coasdinc  and  the 
offshore  barrier  islands  near  Biloxi,  MS,  as  seen  in  Figure  63d.  As  CamiUe 
moves  onto  land  (Figures  63e  and  63f),  the  error  in  computed  storm  surge 
shifts  from  overprediction  to  underprediction  as  winds  come  off  the  land.  In 
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Figure  63.  Contours  of  the  absolute  overprediction  (reds)  and  underprediction  (blues)  errors  in 
the  storm  surge  computed  over  grid  SG01  (Sheet  1  of  3) 
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Figure  63.  (Sheet  2  of  3) 
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Figure  64.  Enlargement  of  the  resolution  provided  by  grid  SG01  in  the  landfall  region  of  Hurri¬ 
cane  Camille 

general,  the  spatial  distribution  of  the  absolute  error  indicates  that  increased 
resolution  around  offshore  islands  and  in  coastal  areas  near  the  landfall  point 
of  Hurricane  Camille  may  reduce  the  extreme  error  magnitudes  observed  in  the 
storm  surge  prediction. 

A  sampling  of  13  storm  surge  hydrographs  taken  at  Alligator  Bayou, 
Pensacola,  Biloxi,  Bay  St.  Louis,  Cat  Island,  Stations  C.6,  C.9,  and  C.13  along 
the  Gulf  of  Mexico  coast.  Stations  S.8  and  S.9  on  the  continental  shelf.  Sta¬ 
tions  SB. 7  and  SB. 8  along  the  shelf  break,  and  at  the  GOM  Pelagic  in  the 
deep  Gulf  of  Mexico  are  shown  in  Figures  65  -  77.  Each  hydrograph  com¬ 
pares  the  storm  surge  elevations  computed  over  grid  SGOl  to  those  computed 
over  grid  CGOl.  Generally,  the  peak  storm  surge  is  overpredicted  with  errors 
ranging  from  5  cm  at  Pensacola,  FL,  to  70  cm  at  Biloxi,  MS.  On  the  right- 
hand  side  of  the  hurricane,  overprediction  of  the  storm  surge  increases  as  sta¬ 
tion  location  nears  the  proximity  of  the  landfall  location  of  Hurricane  Camille. 
The  hydrograph  at  Station  C.9  exhibits  drying  in  the  form  of  negative  surges 
which  are  associated  with  the  left-hand  side  of  the  hurricane.  At  stations  on 
the  continental  shelf  and  at  the  shelf  break,  seen  in  Figures  73  -  76,  the  pealc 
surge  is  represented  with  considerably  more  accuracy.  Overprediction  enors 
are  commonly  less  than  5  cm,  with  a  maximum  overprediction  of  20  m  seen  at 
Station  S.8  in  Figure  73.  Representation  of  the  surge  forerunner  over 
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Figure  65.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Alligator  Bayou,  FL 


Figure  66.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Pensacola,  FL 
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Figure  67.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Biloxi,  MS 
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Figure  68.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Bay  St.  Louis,  MS 
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Figure  69.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Cat  Island,  MS 


Figure  70.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Station  C.6 
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HURRICANE  CAMILLE 


Figure  71 .  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Station  C.9 


Figure  72.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Station  C.13 
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Figure  73.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Station  S.8 


Figure  74.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Station  S.9 
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HURRICANE  CAMILLE 


Figure  75.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Station  SB.7 


Figure  76.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Station  SB.8 
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Figure  77.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  the  GOM  (Gulf  of  Mexico)  Pelagic  station 

grid  SGOl  differs  by  less  than  5  m  with  storm  surge  elevations  computed  over 
grid  CGOl. 

As  predicted  by  the  study  of  grid  structure  in  Chapter  5,  the  most 
significant  errors  are  computed  in  regions  where  grid  refinement  is  less  and  the 
detail  of  the  shoreline  is  more  complex.  Overall,  errors  in  the  storm  surge 
elevations  computed  over  grid  SGOl  are  uniform  and  minimized  throughout 
the  domain.  Further  refinement  in  coastal  areas  near  Biloxi,  MS,  may  improve 
representation  of  the  peak  surge  as  the  hurricane  makes  landfall. 


Hurricane  Kate 

Locations  and  magnitudes  of  the  maximum  storm  surge  prediction  errors 
for  Hurricane  Kate  computed  over  grid  SGOl  are  presented  in  Figure  78. 
Similar  to  the  analysis  for  Hurricane  Camille  forcing,  a  “truth”  solution  is 
defined  as  surge  elevations  computed  over  grid  CGOl.  In  Figure  78,  the  maxi¬ 
mum  absolute  errors  remain  uniform  throughout  the  simulation  and  range 
between  20  and  25  cm  of  overprediction  and  25  and  30  cm  of  underprediction. 
The  largest  relative  errors  in  Figure  78  are  computed  during  the  beginning  of 
the  simulation,  while  Hurricane  Kate  is  in  the  deep  ocean  and  are  a  result  of 
an  under-resolution  of  the  spatial  scale  of  Hurricane  Kate  over  deep  Atlantic 
waters.  At  108  hr,  a  peak  of  the  maximum  relative  error  seen  in  Figure  78 
coincides  with  resonant  mode  generation  prior  to  the  time  of  peak  surge. 
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Figure  78.  Maximum  over-  and  underprediction  errors  in  the  storm  surge 
generated  by  Hurricane  Kate  as  computed  over  grid  SG01 

Spatial  distributions  of  the  absolute  error  associated  with  the  storm  surge 
calculated  over  grid  SGOl  and  compared  to  computations  over  grid  CGOl  are 
presented  at  different  times  in  Figures  79a  -  79f.  In  Figure  79,  the  overpredic¬ 
tion  of  storm  surge  is  shown  in  shades  of  red  and  the  underprediction  in 
shades  of  blue.  As  Hurricane  Kate  moves  into  the  northeast  Gulf  of  Mexico, 
seen  in  Figures  79a  and  79b,  overprediction  errors  along  the  right-  hand  side 
of  the  hurricane  are  evident  in  regions  where  the  discretization  is  relatively 
coarse.  Under-resolved  areas  are  evident  upon  examination  of  the  grid  discret¬ 
ization  in  the  landfall  region  of  Hurricane  Kate  shown  in  Figure  80.  The 
extreme  errors  seen  in  Figure  79c  occur  near  the  islands  offshore  from 
Apalachicola,  FL.  At  150  hr,  approximately  the  time  of  landfall  of  Hurricane 
Kate,  errors  in  the  predicted  storm  surge  near  the  landfall  location  at  Panama 
City,  FL,  shown  in  Figure  79d,  are  noticeable  and  localized  errors  around 
Apalachicola  persist.  Elsewhere  the  storm  surge  is  underpredicted,  particularly 
over  the  less  resolved  portions  of  the  Gulf  of  Mexico  away  from  the  shoreline. 

Following  the  time  of  peak  surge,  underprediction  of  the  water  elevation 
occurs  in  the  northeast  comer  of  the  Gulf  of  Mexico  (Figure  6.25e)  where 
resolution  is  particularly  coarse  relative  to  the  surrounding  coastline.  As  the 
storm  moves  further  inland  and  dissipates  (Figure  6.25f),  errors  in  the  pre¬ 
dicted  Storm  surge  return  to  levels  of  1  to  2  cm  (0.4-.8  in.).  In  general,  errors 
in  the  storm  surge  predicted  for  Hurricane  Kate  forcing  are  quite  low  (i.e.,  less 
than  5  cm  (2.0  in.)).  Only  errors  in  a  few  localized  areas  remain  problematic. 
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Contours  of  the  absolute  overprediction  (reds)  and  underprediction  (blues)  errors  in 
the  storm  surge  computed  over  grid  SG01  (Sheet  1  of  3) 


Figure  79.  (Sheet  2  of  3) 


(c)  144  hr 


(d)  150  hr 
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Figure  80.  Enlargement  of  the  resolution  provided  by  grid  SGOl  in  the  landfall  region  of 
Hurricane  Kate 

A  sampling  of  13  storm  surge  hydrographs  at  Alligator  Bayou,  Destin, 
Panama  City,  Apalachicola,  Carrabelle,  Shell  Point,  and  Cedar  Key  along  the 
Gulf  of  Mexico  coast.  Stations  S.l,  S.2,  S.3,  and  the  outer  Florida  Shelf  station 
on  the  continental  shelf,  and  Stations  SB.l  and  SB. 3  along  the  shelf  break  in 
the  Gulf  of  Mexico  are  shown  in  Figures  81  -  93.  Each  hydrograph  compares 
storm  surge  elevations  computed  over  grid  SGOl  to  those  computed  over 
grid  CGOl.  The  peak  storm  surge  at  all  stations  is  well-predicted,  having 
errors  less  than  4  cm  (1.6  in.)  at  all  stations  except  Apalachicola  where  errors 
of  10  cm  (3.9  in.)  are  computed.  A  majority  of  the  error  in  the  computed 
storm  surge  for  Hurricane  Kate  arises  when  representing  the  resonant  modes, 
which  are  excited  as  Hurricane  Kate  enters  and  moves  through  the  Gulf  of 
Mexico.  Note  that  while  errors  in  the  prediction  of  surge  forerunner  persist, 
error  magnitudes  are  generally  less  than  5  cm  (2.0  in.). 

Storm  surge  heights  computed  over  grid  SGOl  due  to  forcing  from  Hurri¬ 
cane  Kate  are  well-represented  when  compared  to  the  storm  surge  calculated 
for  the  over-refined  “truth”  grid,  CGOl.  Errors  in  the  computed  storm  surge 
relative  to  Hurricane  Kate  are  similar  to  errors  in  the  storm  surge  elevations 
predicted  for  Hurricane  Camille  in  that  they  are  confined  to  very  localized 
areas  where  resolution  is  coarse  and  coastline  detail  is  highly  irregular. 


136 


Chapter  6  Application  of  Domain  Size  and  Gridding  Strategy 


Figure  83.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Panama  City,  FL 


Figure  84.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Apalachicola,  FL 
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NOVEMBER,  1985 


Figure  87.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Cedar  Key,  FL 


Figure  88.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Station  S.1 
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Figure  89.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Station  S.2 


Figure  90.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Station  S.3 
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Figure  91 .  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  the  outer  Florida  Shelf  Station 


Figure  92.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Station  SB.1 
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Figure  93.  Storm  surge  hydrograph  comparing  computations  over  grids  SG01 
and  CG01  at  Station  SB.3 


Conclusions 

A  numerical  storm  surge  model  implementing  a  variably  graded  grid  struc¬ 
ture  in  combination  with  a  large  domain  which  extends  from  the  coastal  region 
to  portions  of  the  deep  ocean  (e.g.  grid  SGOl)  reproduces  weU  the  primary 
storm  surge  elevations  generated  by  Hurricane  Kate  and  Hurricane  CamiUe. 
Errors  in  the  predicted  storm  surge  associated  with  each  hurricane  remain 
fairly  uniform  over  the  domain,  which  is  essential  to  eliminating  tlie  influence 
of  the  discretization  on  storm  surge  computations.  Resolution  over  the  deep 
open  ocean  for  forcing  due  to  Hurricane  Kate  could  be  improved.  However, 
along  the  path  of  Hurricane  Kate  at  entrances  to  the  Gulf  of  Mexico,  where  the 
setup  of  resonant  modes  is  critical,  mesh  spacing  is  generally  less  than  one  half 
the  spatial  scale  of  the  hurricane,  meeting  the  guideline  set  forth  to  reduce 
storm  surge  prediction  errors  over  deep  water. 

For  both  Hurricane  Kate  and  Hurricane  Camille,  the  largest  stonn  surge 
prediction  errors  occur  near  the  shore  at  the  landfall  point  of  the  hurricane. 

For  Hurricane  Kate,  errors  are  highly  localized  in  the  vicinity  of  the  complex 
morphology  associated  with  coastal  inlets  and  islands.  Errors  in  the  predicted 
storm  surge  elevation  relative  to  Hurricane  Camille  forcing  are  significantly 
larger  in  magnitude  and  extend  over  greater  portions  of  the  nearshore  region 
than  errors  produced  by  Hurricane  Kate  forcing.  However,  the  pattern  of  the 
prediction  errors  are  similar  to  those  seen  relative  to  Hurricane  Kate  computa¬ 
tions.  Storm  surge  prediction  errors  associated  with  Hurricane  Camille  are 
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concentrated  around  the  irregular  shorelines  of  barrier  islands  off  the  Missis¬ 
sippi  and  Florida  coasts  and  in  and  around  coastal  inlets.  The  increased  speed 
of  Hurricane  Camille  relative  to  Hurricane  Kate  may  in  fact  amplify  the  errors 
in  under-resolved  regions. 

The  surge  forerunner  generated  within  the  Gulf  of  Mexico  is  not  well  repre¬ 
sented  by  the  variably  graded  discretization  associated  with  grid  SGOl. 
Meteorological  forcing  from  Hurricane  Kate  produces  noticeable  oscillations  of 
the  storm  surge  elevation.  The  discretization  provided  by  grid  SGOl  leads  to 
5-cm  errors  in  the  prediction  of  these  resonant  modes.  The  work  of  Westerink, 
Luettich,  and  Muccino  (1994)  and  Luettich  and  Westerink  (1994)  has  shown 
that  accurate  prediction  of  periodic  forcing,  such  as  that  caused  by  tides, 
requires  additional  resolution  over  the  continental  shelf  break.  The  discretiza¬ 
tion  provided  by  grid  SGOl  concentrates  refinement  at  the  coastline,  so 
periodic  phenomena  requiring  substantial  resolution  at  the  shelf  break  would 
not  be  accurately  captured. 

The  rectangular  basins  used  for  the  grid  structure  study  in  Chapter  5  were 
not  constructed  to  have  resonant  characteristics.  Consequently,  the  need  for 
resolution  at  the  shelf  break  is  not  evident  from  the  error  analyses  of  storm 
surge  predictions  over  the  rectangular  grids.  Significant  refinement  at  the 
coastline  serves  to  accurately  predict  the  primary  surge  but  not  surge  fore¬ 
runner  effects. 

In  contrast  to  forcing  by  Hurricane  Kate,  the  storm  surge  hydrographs  asso¬ 
ciated  with  Hurricane  Camille  indicate  that  the  generation  of  resonant  modes  is 
not  as  prevalent.  This  may  be  caused  in  part  by  the  fact  that  Hurricane 
Camille  enters  the  Gulf  of  Mexico,  at  least  partially,  from  the  Island  of  Cuba 
causing  perhaps  an  incomplete  excitation  of  resonant  modes  in  the  Yucatan 
Channel.  Secondly,  the  rapid  speed  of  Hurricane  Camille  may  not  allow 
proper  setup  of  the  resonant  modes  within  the  Gulf  of  Mexico. 

One  further  point  is  that  the  grid  convergence  techniques  used  to  evaluate 
errors  in  the  computed  storm  surge  relative  to  the  grid  discretization  are  easily 
implemented  and  readily  identify  regions  within  the  grid  discretization  which 
may  require  additional  refinement.  The  use  of  storm  surge  elevations  com¬ 
puted  over  an  over-refined  “truth”  grid  which  has  a  fourfold  refinement  of  the 
study  grid  resolution  is  a  convenient  and  tractable  way  of  assessing  errors  in 
the  predicted  storm  surge  associated  with  a  particular  grid  discretization.  The 
error  analysis  procedure  presented  offers  an  efficient  means  of  evaluating 
storm  surge  prediction  errors  and  consequently  assists  in  the  construction  of  a 
grid  discretization  which  is  both  computationally  efficient  and  yields  accurate 
predictions  of  hurricane  storm  surge  in  coastal  regions. 
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7  Concludinq  Remarks 


Numerical  models  of  continental  margin  waters  are  a  powerful  tool  when 
properly  applied.  Most  coastal  ocean  problems  of  interest  involve  highly  non¬ 
linear  hydrodynamic  processes  which  occur  over  a  wide  range  of  scales.  An 
accurate  model  must  represent  the  complex  coastal  geometries  along  the  shore¬ 
line  and  in  and  around  inlets  and  embayments  as  well  as  incorporate  the 
detailed  bathymetric  variations  of  the  continental  shelf  and  slope  region.  Real¬ 
istic  coastal  ocean  problems  preclude  solution  by  analytical  means.  As  a 
result,  numerical  model  computations  assume  a  primary  role  in  the  study  of  the 
physics  of  continental  margin  waters. 

The  success  of  a  numerical  model  depends  on  the  fonnulation  of  each  of  its 
components  including  the  governing  equations,  numerics,  boundary  conditions, 
grid  discretization,  and  the  computational  domain.  Convergence  studies  are 
undertaken  to  assess  the  perfoniiance  of  model  components  in  simulating  coas¬ 
tal  ocean  processes.  Rigorous  convergence  testing  leads  to  formulations  of  the 
model  components  which  have  a  minimal  impact  on  the  computed  solution. 
Furthermore,  any  remaining  errors  due  to  the  discrete  model  fonn  are  quanti¬ 
fied.  Upon  application  of  a  numerical  model  which  has  been  subjected  to 
extensive  convergence  testing,  computations  can  be  interpreted  with  respect  to 
the  physical  phenomena  simulated  and  any  numerical  artifacts  which  may 
reside  in  the  computed  solution.  Convergence  properties  of  a  numerical  model 
must  be  well  understood  to  obtain  meaningful  results  from  numerical  model 
computations. 

For  the  study  of  stonn  surge  generation  and  propagation  in  the  coastal 
ocean  due  to  hurricane  forcing,  a  numerical  stonn  surge  model  is  implemented. 
Convergence  of  numerical  stonn  surge  models  with  respect  to  the  governing 
equations,  numerics,  and  boundary  condition  fonnulations  is  well  established. 
However,  no  previous  efforts  have  been  made  to  assess  the  effects  of  domain 
size  and  grid  structure  on  numerical  stonn  surge  predictions.  The  work  per- 
fonned  herein  has  quantified  the  influences  of  domain  specification  and  grid 
discretization  on  the  computation  of  stonn  surge  elevation.  As  a  consequence, 
a  series  of  guidelines  is  set  forth  to  direct  the  construction  of  the  computational 
domain  and  its  discretization  so  that  accurate  stonn  surge  predictions  result. 

Development  of  a  hydrodynamic  model  has  not  been  the  focus  of  the  work 
presented,  but  the  success  of  the  stomi  surge  convergence  studies  undertaken 
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does  depend  on  selection  of  an  accurate,  efficient,  and  flexible  numerical  storm 
surge  model.  A  finite  element,  two-dimensional  hydrodynamic  model, 
ADCIRC-2DDI,  based  on  the  generalized  wave  continuity  formulation  of  the 
depth-averaged  conservation  equations  fulfills  the  accuracy,  efficiency,  and 
grid  flexibility  requirements  set  forth  and  has  provided  an  appropriate  frame¬ 
work  for  modeling  the  storm  surge  response  of  the  coastal  ocean.  Conver¬ 
gence  of  this  model  has  already  been  proven  with  regard  to  the  governing 
equations,  numerics,  and  boundary  condition  formulation. 

The  domain  size  sensitivity  study  conducted  establishes  a  relationship 
between  domain  size,  boundary  condition  specification,  and  the  resulting 
physics  associated  with  hurricane  storm  surge  generation  in  continental  margin 
waters.  Comparisons  of  stonn  surge  computations  over  three  domain  sizes 
subject  to  two  different  open  ocean  boundary  forcings  form  the  core  of  this 
study. 

Results  clearly  illustrate  that  a  large  stonn  surge  model  domain,  which 
includes  the  western  North  Atlantic  Ocean,  the  Gulf  of  Mexico,  and  the  Carib¬ 
bean  Sea,  leads  to  convergent  predictions  of  both  the  primary  stonn  surge  and 
surge  forerunner  effects.  The  inclusion  of  contiguous  basins  within  the  large 
domain  allows  proper  setup  of  basin  resonant  modes  and  facilitates  the  realistic 
propagation  of  storm  surge  onto  the  continental  shelf  where  development  of 
stonn  surge  is  most  critical.  The  main  advantage  of  the  large  domain  is  that 
open  boundaries  lie  within  the  deep  Atlantic  Ocean  and  are  far  from  the  intri¬ 
cate  processes  occurring  on  and  around  the  continental  shelf 

Study  of  domain  size  demonstrates  that  the  commonly  used  computational 
domain,  one  that  is  situated  on  the  continental  shelf  and  whose  size  is  limited 
relative  to  the  size  of  the  storm,  significantly  underestimates  the  primary  storm 
surge  response.  The  large  expanse  of  cross-shelf  open  ocean  boundaries, 
located  in  the  vicinity  of  significant  storm  surge  generation,  prevents  accurate 
water  elevation  specifications  at  the  open  boundaries  of  this  small  domain. 
Furthermore,  a  domain  encompassing  a  resonant  basin,  such  as  the  Gulf  of 
Mexico,  which  is  quite  sensitive  to  boundary  forcing  functions,  may  not  cap¬ 
ture  the  physics  associated  with  the  hurricane  forerunner. 

Having  established  the  convergence  of  the  stonn  surge  model  with  respect 
to  domain  size,  an  investigation  into  the  influence  of  the  domain  discretization 
on  storm  surge  computations  proceeds.  Comparisons  of  stonn  surge  elevations 
computed  over  14  grids  subject  to  4  synthetic  hurricane  forcings  substantiates 
the  need  for  significant  resolution  in  near-coastal  regions  and  of  the  shoreline 
itself  For  the  hurricanes  considered,  storm  surge  generation  is  shown  to  inten¬ 
sify  for  storms  of  large  spatial  scale,  those  moving  rapidly,  and  for  storms 
which  make  landfall  at  the  coast.  Accurate  prediction  of  the  stonn  surge  for 
hurricanes  having  these  characteristics  requires  even  higher  levels  of  grid 
refinement  in  the  near-shore  region.  Complex  coastline  geometries  also 
increase  stonn  surge  generation  and  can  lead  to  substantial  error  in  predicted 
storm  surge  elevations  when  shoreline  detail  is  misrepresented.  In  the  deep 
ocean,  the  inverted  barometer  effect  is  captured  by  grid  refinements  of 
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approximately  one  half  the  spatial  scale  of  the  hurricane  as  measured  by  the 
radius  to  maximum  wind.  The  disparity  in  discretization  requirements  over  a 
large  domain  render  the  variably  graded,  unstructured  discretization  as  optimal. 

The  variably  graded  grid  structure  yields  low  uniform  errors  in  the  storm  surge 
prediction  throughout  the  domain  and  minimizes  computational  effort. 

Application  of  the  guidelines  for  domain  size  and  grid  discretization  strat¬ 
egy  to  two  historical  hurricanes  (Hurricane  Camille  (1969)  and  Hurricane  Kate 
(1985))  reinforces  the  findings  that  high  levels  of  resolution  in  the  nearshore 
region  and  in  areas  of  complex  coastal  geometry  are  important  for  accurate 
prediction  of  primary  storm  surge.  This  is  particularly  noticeable  for  storm 
surge  generation  associated  with  the  intense  hurricane  forcing  of  Hurricane 
CamiUe.  These  applications  indicate  that  the  resolution  provided  at  the  coast¬ 
line  by  the  study  grid  is  not  sufficient  for  proper  representation  of  hurricane 
forerunner/resonant  mode  behavior  and  that  additional  resolution  over  the 
continental  shelf  break  may  be  necessary. 

In  the  context  of  these  applications,  a  procedure  to  assess  the  performance 
of  a  particular  grid  discretization  is  demonstrated.  The  analysis  techniques 
implemented  are  easily  applied  and  extremely  effective  in  identifying  portions 
of  the  grid  discretization  which  may  be  under-resolved.  The  convergence 
analysis  presented  through  these  applications  demonstrates  the  form  of  testing 
which  should  be  conducted  prior  to  any  numerical  stomi  surge  modeling  effort. 

In  summary,  recall  that  hurricanes  are  large-scale  phenomena  whose  winds 
cover  significant  areal  regions.  As  such,  a  domain  which  has  an  areal  extent 
much  greater  than  the  scale  of  the  hurricane  is  likely  to  more  realistically 
capture  the  physics  associated  with  stonn  surge  generation  and  propagation 
without  requiring  a  detailed  a  priori  knowledge  of  the  hydrodynamics  at  the 
open  ocean  boundaries  of  the  domain.  Thus,  the  complex  oceanic  response 
generated  by  hurricane  wind  forcing  is  best  represented  when  modeled  near  its 
inception  in  the  deep  Atlantic  ocean  and  tracked  naturally  through  contiguous 
basins,  onto  the  continental  shelf,  and  into  shallow  coastal  regions.  Further¬ 
more,  the  range  of  grid  spacings  required  to  accurately  represent  the  stonn 
surge  throughout  the  domain  is  best  accommodated  through  the  use  of  a  vari¬ 
ably  graded  grid  structure.  The  graded  grid  discretization  has  significant 
refinement  at  the  coastline  and  in  nearshore  regions  and  resolution  to  one  half 
the  spatial  scale  of  the  storm  in  the  deep  ocean.  Using  these  guidelines  for 
selection  of  a  domain  size  and  construction  of  a  grid  discretization,  errors  in 
the  prediction  of  stonn  surge  over  the  domain  remain  low  and  uniform. 

Finally,  flexibility  of  the  finite  element  method  ensures  that  storm  surge  can  be 
computed  efficiently  over  large  domains  with  grid  spacings  ranging  over  three 
to  four  orders  of  magnitude. 

Aside  from  the  reported  findings  of  the  stonn  surge  convergence  studies, 
several  issues  arise  which  require  further  investigation.  Accurate  representa¬ 
tion  of  resonant  modes  generated  by  hurricane  forcing  is  beyond  the  scope  of 
the  work  presented  but  remains  an  important  area  of  study.  The  work  of 
Westerink,  Luettich,  and  Muccino  (1994)  and  Luettich  and  Westerink  (in 
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preparation)  demonstrates  that  additional  resolution  over  the  shelf  break  is 
necessary  to  properly  represent  tidal  forcing.  Refinement  of  the  shelf  break 
may  also  improve  the  accurate  propagation  of  resonant  modes  onto  the  conti¬ 
nental  shelf.  Furthermore,  storm  surge  simulations  of  Hurricane  Kate  and 
Hurricane  Camille  have  shown  vastly  different  resonant  characteristics,  indicat¬ 
ing  a  complex  relationship  between  the  influence  of  the  basin  and  the  hurri¬ 
cane  forcing  on  the  excitation  of  resonant  modes.  Resonant  modes/surge  fore¬ 
runner  are  an  important  component  of  the  overall  storm  surge  elevation  but 
their  generation  is  not  well  understood. 

AU  of  the  convergence  studies  presented  are  conducted  relative  to  the  pre¬ 
diction  of  water  elevations.  Clearly,  analogous  studies  are  also  needed  to 
examine  the  effects  of  domain  size  and,  in  particular,  grid  structure  on  velocity 
computations.  Moreover,  all  computations  were  performed  using  two- 
dimensional,  depth-averaged  equations.  Three-dimensional  storm  surge  simula¬ 
tions  are  a  logical  progression  to  determine  under  what  conditions  use  of  the 
depth-averaged  equations  in  storm  surge  prediction  is  inappropriate. 

Finally,  use  of  the  HURWIN  wind  model  for  the  generation  of  hurricane 
wind  and  pressure  fields  reinforces  the  sensitivity  of  the  specification  of  mete¬ 
orological  forcing  on  the  hydrodynamic  model  response.  Sensitivity  and  veri¬ 
fication  studies  undertaken  with  respect  to  the  meteorological  forcing  are  in 
order,  as  weU  as  a  quantification  of  resolution  requirements  for  the  wind  forc¬ 
ing  in  relation  to  the  hydrodynamic  model  grid  discretization. 
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